/* * Copyright © 2003-2022 Dynare Team * * This file is part of Dynare. * * Dynare is free software: you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by * the Free Software Foundation, either version 3 of the License, or * (at your option) any later version. * * Dynare is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * GNU General Public License for more details. * * You should have received a copy of the GNU General Public License * along with Dynare. If not, see . */ #ifndef _MODELTREE_HH #define _MODELTREE_HH #include #include #include #include #include #include #include #include #include #include "DataTree.hh" #include "EquationTags.hh" #include "ExtendedPreprocessorTypes.hh" #include "Bytecode.hh" using namespace std; // Helper to convert a vector into a tuple template auto vectorToTupleHelper(const vector &v, index_sequence) { return tuple{v[Indices]...}; } template auto vectorToTuple(const vector &v) { assert(v.size() >= N); return vectorToTupleHelper(v, make_index_sequence()); } //! Vector describing equations: BlockSimulationType, if BlockSimulationType == EVALUATE_s then a expr_t on the new normalized equation using equation_type_and_normalized_equation_t = vector>; //! Vector describing variables: max_lag in the block, max_lead in the block using lag_lead_vector_t = vector>; //! Shared code for static and dynamic models class ModelTree : public DataTree { friend class DynamicModel; friend class StaticModel; public: // Set via the `compiler` command string user_set_add_flags, user_set_subst_flags, user_set_add_libs, user_set_subst_libs, user_set_compiler; protected: /* * ************** BEGIN ************** * The following structures keep track of the model equations and must all be updated * when adding or removing an equation. Hence, if a new parallel structure is added * in the future, it must be maintained whereever these structures are updated * See in particular methods clearEquations(), replaceMyEquations() and * computeRamseyPolicyFOCs() of DynamicModel class. * NB: This message added with the introduction of the `exclude_eqs` option, hence * that's a place to update future structures. */ //! Stores declared and generated auxiliary equations vector equations; /* Stores line numbers of declared equations; undefined in some cases (e.g. auxiliary equations) */ vector> equations_lineno; //! Stores equation tags EquationTags equation_tags; /* * ************** END ************** */ //! Only stores generated auxiliary equations, in an order meaningful for evaluation /*! These equations only contain the definition of auxiliary variables, and may diverge from those in the main model (equations), if other model transformations applied subsequently. This is not a problem, since aux_equations is only used for regenerating the values of auxiliaries given the others. For example, such a divergence appears when there is an expectation operator in a ramsey model, see tests/optimal_policy/nk_ramsey_expectation.mod */ vector aux_equations; //! Maximum order at which (endogenous) derivatives have been computed int computed_derivs_order{0}; //! Stores derivatives /*! Index 0 is not used, index 1 contains first derivatives, ... For each derivation order, stores a map whose key is a vector of integer: the first integer is the equation index, the remaining ones are the derivation IDs of variables (in non-decreasing order, to avoid storing symmetric elements several times). Only non-zero derivatives are stored. */ vector, expr_t>> derivatives; //! Number of non-zero derivatives /*! Index 0 is not used, index 1 contains number of non-zero first derivatives, ... */ vector NNZDerivatives; //! Derivatives with respect to parameters /*! The key of the outer map is a pair (derivation order w.r.t. endogenous, derivation order w.r.t. parameters). For e.g., { 1, 2 } corresponds to the jacobian differentiated twice w.r.t. to parameters. In inner maps, the vector of integers consists of: the equation index, then the derivation IDs of endogenous (in non-decreasing order), then the IDs of parameters (in non-decreasing order)*/ map, map, expr_t>> params_derivatives; //! Used model local variables, that will be treated as temporary terms /*! See the comments in ModelTree::computeTemporaryTerms() */ map temporary_terms_mlv; //! Temporary terms for residuals and derivatives /*! Index 0 is temp. terms of residuals, index 1 for first derivatives, ... */ vector temporary_terms_derivatives; //! Stores, for each temporary term, its index in the MATLAB/Julia vector temporary_terms_idxs_t temporary_terms_idxs; //! Temporary terms for parameter derivatives, under a disaggregated form /*! The pair of integers is to be interpreted as in param_derivatives */ map, temporary_terms_t> params_derivs_temporary_terms; //! Stores, for each temporary term in param. derivs, its index in the MATLAB/Julia vector temporary_terms_idxs_t params_derivs_temporary_terms_idxs; //! Trend variables and their growth factors map trend_symbols_map; //! for all trends; the boolean is true if this is a log-trend, false otherwise using nonstationary_symbols_map_t = map>; //! Nonstationary variables and their deflators nonstationary_symbols_map_t nonstationary_symbols_map; /* Maps indices of equations in the block-decomposition order into original equation IDs */ vector eq_idx_block2orig; /* Maps indices of (endogenous) variables in the block-decomposition order into original type-specific endogenous IDs */ vector endo_idx_block2orig; /* Maps original variable and equation indices into the block-decomposition order. Set by updateReverseVariableEquationOrderings() */ vector eq_idx_orig2block, endo_idx_orig2block; //! Vector describing equations: BlockSimulationType, if BlockSimulationType == EVALUATE_s then a expr_t on the new normalized equation equation_type_and_normalized_equation_t equation_type_and_normalized_equation; /* Stores derivatives of each block w.r.t. endogenous that belong to it. The tuple is: equation number (inside the block), variable number (inside the block), lead/lag */ vector, expr_t>> blocks_derivatives; class BlockInfo { public: BlockSimulationType simulation_type; int first_equation; // Stores a block-ordered equation ID int size{0}; int mfs_size{0}; // Size of the minimal feedback set bool linear{true}; // Whether the block is linear in endogenous variable int n_static{0}, n_forward{0}, n_backward{0}, n_mixed{0}; int max_endo_lag{0}, max_endo_lead{0}; // Maximum lag/lead on endos that appear in and *that belong to* the block int max_other_endo_lag{0}, max_other_endo_lead{0}; // Maximum lag/lead on endos that appear in but do not belong to the block int max_exo_lag{0}, max_exo_lead{0}; int max_exo_det_lag{0}, max_exo_det_lead{0}; int max_lag{0}, max_lead{0}; // The max over all endo/exo variables int getRecursiveSize() const { return size - mfs_size; }; }; // Stores various informations on the blocks vector blocks; // Maps endogenous type-specific IDs to the block number to which it belongs vector endo2block; /* Maps (original) equation number to the block number to which it belongs. It verifies: ∀i, eq2block[endo2eq[i]] = endo2block[i] */ vector eq2block; /* Temporary terms for block decomposed models. - the outer vector has as many elements as there are blocks in the model - the inner vector has as many elements as there are equations in the block, plus a last one which contains the temporary terms for derivatives It’s necessary to track temporary terms per equation, because some equations are evaluated instead of solved, and an equation E1 may depend on the value of an endogenous Y computed by a previously evaluated equation E2; in this case, if some temporary term TT of equation E2 contains Y, then TT needs to be computed after E1, but before E2. */ vector> blocks_temporary_terms; /* Stores, for each temporary term in block decomposed models, its index in the vector of all temporary terms */ temporary_terms_idxs_t blocks_temporary_terms_idxs; //! Computes derivatives /*! \param order the derivation order \param vars the derivation IDs w.r.t. which compute the derivatives */ void computeDerivatives(int order, const set &vars); //! Computes derivatives of the Jacobian and Hessian w.r. to parameters void computeParamsDerivatives(int paramsDerivsOrder); //! Computes temporary terms (for all equations and derivatives) void computeTemporaryTerms(bool is_matlab, bool no_tmp_terms); //! Computes temporary terms per block void computeBlockTemporaryTerms(); /* Add additional temporary terms for a given block. This method is called by computeBlockTemporaryTerms(). It does nothing by default, but is meant to be overriden by subclasses (actually by DynamicModel, who needs extra temporary terms for derivatives w.r.t. exogenous and other endogenous) */ virtual void additionalBlockTemporaryTerms(int blk, vector> &blocks_temporary_terms, map> &reference_count) const; //! Computes temporary terms for the file containing parameters derivatives void computeParamsDerivativesTemporaryTerms(); //! Writes temporary terms template void writeTemporaryTerms(const temporary_terms_t &tt, temporary_terms_t &temp_term_union, const temporary_terms_idxs_t &tt_idxs, ostream &output, deriv_node_temp_terms_t &tef_terms) const; void writeJsonTemporaryTerms(const temporary_terms_t &tt, temporary_terms_t &temp_term_union, ostream &output, deriv_node_temp_terms_t &tef_terms, const string &concat) const; //! Writes temporary terms in bytecode template void writeBytecodeTemporaryTerms(const temporary_terms_t &tt, temporary_terms_t &temporary_terms_union, BytecodeWriter &code_file, deriv_node_temp_terms_t &tef_terms) const; /* Adds information for (non-block) bytecode simulation in a separate .bin file. Returns the number of first derivatives w.r.t. endogenous variables */ int writeBytecodeBinFile(const string &filename, bool is_two_boundaries) const; //! Adds per-block information for bytecode simulation in a separate .bin file int writeBlockBytecodeBinFile(ofstream &bin_file, int block) const; //! Fixes output when there are more than 32 nested parens, Issue #1201 void fixNestedParenthesis(ostringstream &output, map &tmp_paren_vars, bool &message_printed) const; //! Tests if string contains more than 32 nested parens, Issue #1201 bool testNestedParenthesis(const string &str) const; template void writeModelLocalVariableTemporaryTerms(temporary_terms_t &temp_term_union, const temporary_terms_idxs_t &tt_idxs, ostream &output, deriv_node_temp_terms_t &tef_terms) const; //! Writes model equations template void writeModelEquations(ostream &output, const temporary_terms_t &temporary_terms) const; // Returns outputs for derivatives and temporary terms at each derivation order template pair, vector> writeModelFileHelper() const; // Writes per-block residuals and temporary terms (incl. for derivatives) template void writePerBlockHelper(int blk, ostream &output, temporary_terms_t &temporary_terms) const; /* Helper for writing derivatives w.r.t. parameters. Returns { tt, rp, gp, rpp, gpp, hp, g3p }. g3p is empty if requesting a static output type. */ template tuple writeParamsDerivativesFileHelper() const; // Helper for writing bytecode (without block decomposition) template void writeBytecodeHelper(BytecodeWriter &code_file) const; // Helper for writing blocks in bytecode template void writeBlockBytecodeHelper(BytecodeWriter &code_file, int block) const; /* Write additional derivatives w.r.t. to exogenous, exogenous det and other endo in block+bytecode mode. Does nothing by default, but overriden by DynamicModel which needs those. */ virtual void writeBlockBytecodeAdditionalDerivatives(BytecodeWriter &code_file, int block, const temporary_terms_t &temporary_terms_union, const deriv_node_temp_terms_t &tef_terms) const; /* Helper for writing JSON output for residuals and derivatives. Returns mlv and derivatives output at each derivation order. */ template pair> writeJsonComputingPassOutputHelper(bool writeDetails) const; /* Helper for writing JSON derivatives w.r.t. parameters. Returns { mlv, tt, rp, gp, rpp, gpp, hp, g3p }. g3p is empty if requesting a static output type. */ template tuple writeJsonParamsDerivativesHelper(bool writeDetails) const; //! Writes JSON model equations //! if residuals = true, we are writing the dynamic/static model. //! Otherwise, just the model equations (with line numbers, no tmp terms) void writeJsonModelEquations(ostream &output, bool residuals) const; /* Writes JSON model local variables. Optionally put the external function variable calls into TEF terms */ void writeJsonModelLocalVariables(ostream &output, bool write_tef_terms, deriv_node_temp_terms_t &tef_terms) const; //! Writes model equations in bytecode template void writeBytecodeModelEquations(BytecodeWriter &code_file, const temporary_terms_t &temporary_terms, const deriv_node_temp_terms_t &tef_terms) const; //! Writes LaTeX model file void writeLatexModelFile(const string &mod_basename, const string &latex_basename, ExprNodeOutputType output_type, bool write_equation_tags) const; //! Sparse matrix of double to store the values of the static Jacobian /*! First index is equation number, second index is endogenous type specific ID */ using jacob_map_t = map, double>; //! Normalization of equations, as computed by computeNonSingularNormalization() /*! Maps endogenous type specific IDs to equation numbers */ vector endo2eq; /* Compute a pseudo-Jacobian whose all elements are either zero or one, depending on whether the variable symbolically appears in the equation */ jacob_map_t computeSymbolicJacobian() const; // Compute {var,eq}_idx_orig2block from {var,eq}_idx_block2orig void updateReverseVariableEquationOrderings(); //! Compute the matching between endogenous and variable using the jacobian contemporaneous_jacobian /*! \param contemporaneous_jacobian Jacobian used as an incidence matrix: all elements declared in the map (even if they are zero), are used as vertices of the incidence matrix \return True if a complete normalization has been achieved */ bool computeNormalization(const jacob_map_t &contemporaneous_jacobian, bool verbose); //! Try to compute the matching between endogenous and variable using a decreasing cutoff /*! Applied to the jacobian contemporaneous_jacobian and stop when a matching is found. If no matching is found using a strictly positive cutoff, then a zero cutoff is applied (i.e. use a symbolic normalization); in that case, the method adds zeros in the jacobian matrices to reflect all the edges in the symbolic incidence matrix. If no matching is found with a zero cutoff, an error message is printed. The resulting normalization is stored in endo2eq. Returns a boolean indicating success. */ bool computeNonSingularNormalization(const jacob_map_t &contemporaneous_jacobian); //! Evaluate the jacobian (w.r.t. endogenous) and suppress all the elements below the cutoff /*! Returns the contemporaneous_jacobian. Elements below the cutoff are discarded. External functions are evaluated to 1. */ jacob_map_t evaluateAndReduceJacobian(const eval_context_t &eval_context) const; /* Search the equations and variables belonging to the prologue and the epilogue of the model. Initializes “eq_idx_block2orig” and “endo_idx_block2orig”. Returns the sizes of the prologue and epilogue. */ pair computePrologueAndEpilogue(); //! Determine the type of each equation of model and try to normalize the unnormalized equation void equationTypeDetermination(const map, expr_t> &first_order_endo_derivatives, int mfs); /* Fills the max lags/leads and n_{static,mixed,forward,backward} fields of a given block. Needs the fields size and first_equation. */ void computeDynamicStructureOfBlock(int blk); /* Fills the simulation_type field of a given block. Needs the fields size, max_endo_lag and max_endo_lead. */ void computeSimulationTypeOfBlock(int blk); /* Compute the block decomposition and for a non-recusive block find the minimum feedback set Initializes the “blocks”, “endo2block” and “eq2block” structures. */ void computeBlockDecomposition(int prologue, int epilogue); /* Reduce the number of block by merging the same type of equations in the prologue and the epilogue */ void reduceBlockDecomposition(); /* The 1st output gives, for each equation (in original order) the (max_lag, max_lead) across all endogenous that appear in the equation and that belong to the same block (i.e. those endogenous are solved in the same block). The 2nd output gives, for each type-specific endo IDs, its (max_lag, max_lead) across all its occurences inside the equations of the block to which it belongs. */ pair getVariableLeadLagByBlock() const; //! Print an abstract of the block structure of the model void printBlockDecomposition() const; //! Determine for each block if it is linear or not void determineLinearBlocks(); //! Return the type of equation belonging to the block EquationType getBlockEquationType(int blk, int eq) const { return equation_type_and_normalized_equation[eq_idx_block2orig[blocks[blk].first_equation+eq]].first; }; //! Return true if the equation has been normalized bool isBlockEquationRenormalized(int blk, int eq) const { return equation_type_and_normalized_equation[eq_idx_block2orig[blocks[blk].first_equation + eq]].first == EquationType::evaluateRenormalized; }; //! Return the expr_t of equation belonging to the block BinaryOpNode * getBlockEquationExpr(int blk, int eq) const { return equations[eq_idx_block2orig[blocks[blk].first_equation + eq]]; }; //! Return the expr_t of renormalized equation belonging to the block BinaryOpNode * getBlockEquationRenormalizedExpr(int blk, int eq) const { return equation_type_and_normalized_equation[eq_idx_block2orig[blocks[blk].first_equation + eq]].second; }; //! Return the original number of equation belonging to the block int getBlockEquationID(int blk, int eq) const { return eq_idx_block2orig[blocks[blk].first_equation + eq]; }; //! Return the original number of variable belonging to the block int getBlockVariableID(int blk, int var) const { return endo_idx_block2orig[blocks[blk].first_equation + var]; }; //! Return the position of an equation (given by its original index) inside its block int getBlockInitialEquationID(int blk, int eq) const { return eq_idx_orig2block[eq] - blocks[blk].first_equation; }; //! Return the position of a variable (given by its original index) inside its block int getBlockInitialVariableID(int blk, int var) const { return endo_idx_orig2block[var] - blocks[blk].first_equation; }; //! Initialize equation_reordered & variable_reordered void initializeVariablesAndEquations(); //! Returns the 1st derivatives w.r.t. endogenous in a different format /*! Returns a map (equation, type-specific ID, lag) → derivative. Assumes that derivatives have already been computed. */ map, expr_t> collectFirstOrderDerivativesEndogenous(); /* Get column number within Jacobian of a given block. “var” is the block-specific endogenous variable index. */ virtual int getBlockJacobianEndoCol(int blk, int var, int lag) const = 0; // Returns a human-readable string describing the model class (e.g. “dynamic model”…) virtual string modelClassName() const = 0; private: //! Internal helper for the copy constructor and assignment operator /*! Copies all the structures that contain ExprNode*, by the converting the pointers into their equivalent in the new tree */ void copyHelper(const ModelTree &m); //! Returns the name of the MATLAB architecture given the extension used for MEX files static string matlab_arch(const string &mexext); #ifdef __APPLE__ //! Finds a suitable GCC compiler on macOS static string findGccOnMacos(const string &mexext); #endif //! Compiles a MEX file void compileMEX(const string &basename, const string &funcname, const string &mexext, const vector &src_files, const filesystem::path &matlabroot, const filesystem::path &dynareroot) const; public: ModelTree(SymbolTable &symbol_table_arg, NumericalConstants &num_constants_arg, ExternalFunctionsTable &external_functions_table_arg, bool is_dynamic_arg = false); protected: ModelTree(const ModelTree &m); ModelTree &operator=(const ModelTree &m); public: //! Absolute value under which a number is considered to be zero double cutoff{1e-15}; //! Compute the minimum feedback set /*! 0 : all endogenous variables are considered as feedback variables 1 : the variables belonging to non normalized equation are considered as feedback variables 2 : the variables belonging to a non linear equation are considered as feedback variables 3 : the variables belonging to a non normalizable non linear equation are considered as feedback variables default value = 0 */ int mfs{0}; //! Declare a node as an equation of the model; also give its line number void addEquation(expr_t eq, optional lineno); //! Declare a node as an equation of the model, also giving its tags void addEquation(expr_t eq, optional lineno, const map &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 int equation_number() const; //! Adds a trend variable with its growth factor void addTrendVariables(const vector &trend_vars, expr_t growth_factor) noexcept(false); //! Adds a nonstationary variables with their (common) deflator void addNonstationaryVariables(const vector &nonstationary_vars, bool log_deflator, expr_t deflator) noexcept(false); //! Is a given variable non-stationary? bool isNonstationary(int symb_id) const; void set_cutoff_to_zero(); /*! Reorder auxiliary variables so that they appear in recursive order in set_auxiliary_variables.m and dynamic_set_auxiliary_series.m */ void reorderAuxiliaryEquations(); //! Find equations of the form “variable=constant”, excluding equations with “mcp” tag (see dynare#1697) void findConstantEquationsWithoutMcpTag(map &subst_table) const; /* Given an expression, searches for the first equation that has exactly this expression on the LHS, and returns the RHS of that equation. If no such equation can be found, throws an ExprNode::MatchFailureExpression */ expr_t getRHSFromLHS(expr_t lhs) const; //! Returns all the equation tags associated to an equation map getEquationTags(int eq) const { return equation_tags.getTagsByEqn(eq); } //! Returns the vector of non-zero derivative counts const vector & getNNZDerivatives() const { return NNZDerivatives; } //! Returns the vector of temporary terms derivatives const vector & getTemporaryTermsDerivatives() const { return temporary_terms_derivatives; } //!Returns the maximum order of computed derivatives int getComputedDerivsOrder() const { return computed_derivs_order; } static string BlockSim(BlockSimulationType type) { switch (type) { case BlockSimulationType::evaluateForward: return "EVALUATE FORWARD "; case BlockSimulationType::evaluateBackward: return "EVALUATE BACKWARD "; case BlockSimulationType::solveForwardSimple: return "SOLVE FORWARD SIMPLE "; case BlockSimulationType::solveBackwardSimple: return "SOLVE BACKWARD SIMPLE "; case BlockSimulationType::solveTwoBoundariesSimple: return "SOLVE TWO BOUNDARIES SIMPLE "; case BlockSimulationType::solveForwardComplete: return "SOLVE FORWARD COMPLETE "; case BlockSimulationType::solveBackwardComplete: return "SOLVE BACKWARD COMPLETE "; case BlockSimulationType::solveTwoBoundariesComplete: return "SOLVE TWO BOUNDARIES COMPLETE"; default: return "UNKNOWN "; } } }; template void ModelTree::writeTemporaryTerms(const temporary_terms_t &tt, temporary_terms_t &temp_term_union, const temporary_terms_idxs_t &tt_idxs, ostream &output, deriv_node_temp_terms_t &tef_terms) const { for (auto it : tt) { if (dynamic_cast(it)) it->writeExternalFunctionOutput(output, output_type, temp_term_union, tt_idxs, tef_terms); it->writeOutput(output, output_type, tt, tt_idxs, tef_terms); output << " = "; it->writeOutput(output, output_type, temp_term_union, tt_idxs, tef_terms); if constexpr(isCOutput(output_type) || isMatlabOutput(output_type)) output << ";"; output << endl; temp_term_union.insert(it); } } template void ModelTree::writeModelLocalVariableTemporaryTerms(temporary_terms_t &temp_term_union, const temporary_terms_idxs_t &tt_idxs, ostream &output, deriv_node_temp_terms_t &tef_terms) const { temporary_terms_t tto; for (const auto &[mlv, value] : temporary_terms_mlv) tto.insert(mlv); for (const auto &[mlv, value] : temporary_terms_mlv) { value->writeExternalFunctionOutput(output, output_type, temp_term_union, tt_idxs, tef_terms); if constexpr(isJuliaOutput(output_type)) output << " const "; mlv->writeOutput(output, output_type, tto, tt_idxs, tef_terms); output << " = "; value->writeOutput(output, output_type, temp_term_union, tt_idxs, tef_terms); if constexpr(isCOutput(output_type) || isMatlabOutput(output_type)) output << ";"; output << endl; /* We put in temp_term_union the VariableNode corresponding to the MLV, not its definition, so that when equations use the MLV, T(XXX) is printed instead of the MLV name */ temp_term_union.insert(mlv); } } template void ModelTree::writeModelEquations(ostream &output, const temporary_terms_t &temporary_terms) const { for (int eq {0}; eq < static_cast(equations.size()); eq++) { BinaryOpNode *eq_node { equations[eq] }; expr_t lhs { eq_node->arg1 }, rhs { eq_node->arg2 }; // Test if the right hand side of the equation is empty. double vrhs {1.0}; try { vrhs = rhs->eval({}); } catch (ExprNode::EvalException &e) { } if (vrhs != 0) // The right hand side of the equation is not empty ==> residual=lhs-rhs; if constexpr(isJuliaOutput(output_type)) { output << " residual" << LEFT_ARRAY_SUBSCRIPT(output_type) << eq + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << " = ("; lhs->writeOutput(output, output_type, temporary_terms, temporary_terms_idxs); output << ") - ("; rhs->writeOutput(output, output_type, temporary_terms, temporary_terms_idxs); output << ")" << endl; } else { output << "lhs = "; lhs->writeOutput(output, output_type, temporary_terms, temporary_terms_idxs); output << ";" << endl << "rhs = "; rhs->writeOutput(output, output_type, temporary_terms, temporary_terms_idxs); output << ";" << endl << "residual" << LEFT_ARRAY_SUBSCRIPT(output_type) << eq + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << " = lhs - rhs;" << endl; } else // The right hand side of the equation is empty ==> residual=lhs; { output << "residual" << LEFT_ARRAY_SUBSCRIPT(output_type) << eq + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << " = "; lhs->writeOutput(output, output_type, temporary_terms, temporary_terms_idxs); output << ";" << endl; } } } template pair, vector> ModelTree::writeModelFileHelper() const { vector d_output(derivatives.size()); // Derivatives output (at all orders, including 0=residual) vector tt_output(derivatives.size()); // Temp terms output (at all orders) deriv_node_temp_terms_t tef_terms; temporary_terms_t temp_term_union; writeModelLocalVariableTemporaryTerms(temp_term_union, temporary_terms_idxs, tt_output[0], tef_terms); writeTemporaryTerms(temporary_terms_derivatives[0], temp_term_union, temporary_terms_idxs, tt_output[0], tef_terms); writeModelEquations(d_output[0], temp_term_union); // Writing Jacobian if (!derivatives[1].empty()) { writeTemporaryTerms(temporary_terms_derivatives[1], temp_term_union, temporary_terms_idxs, tt_output[1], tef_terms); for (const auto &[indices, d1] : derivatives[1]) { auto [eq, var] = vectorToTuple<2>(indices); d_output[1] << "g1" << LEFT_ARRAY_SUBSCRIPT(output_type); if constexpr(isMatlabOutput(output_type) || isJuliaOutput(output_type)) d_output[1] << eq + 1 << "," << getJacobianCol(var) + 1; else d_output[1] << eq + getJacobianCol(var)*equations.size(); d_output[1] << RIGHT_ARRAY_SUBSCRIPT(output_type) << "="; d1->writeOutput(d_output[1], output_type, temp_term_union, temporary_terms_idxs, tef_terms); d_output[1] << ";" << endl; } } // Write derivatives for order ≥ 2 for (size_t i = 2; i < derivatives.size(); i++) if (!derivatives[i].empty()) { writeTemporaryTerms(temporary_terms_derivatives[i], temp_term_union, temporary_terms_idxs, tt_output[i], tef_terms); /* When creating the sparse matrix (in MATLAB or C mode), since storage is in column-major order, output the first column, then the second, then the third. This gives a significant performance boost in use_dll mode (at both compilation and runtime), because it facilitates memory accesses and expression reusage. */ ostringstream i_output, j_output, v_output; for (int k{0}; // Current line index in the 3-column matrix const auto &[vidx, d] : derivatives[i]) { int eq{vidx[0]}; int col_idx{0}; for (size_t j = 1; j < vidx.size(); j++) { col_idx *= getJacobianColsNbr(); col_idx += getJacobianCol(vidx[j]); } if constexpr(isJuliaOutput(output_type)) { d_output[i] << " g" << i << "[" << eq + 1 << "," << col_idx + 1 << "] = "; d->writeOutput(d_output[i], output_type, temp_term_union, temporary_terms_idxs, tef_terms); d_output[i] << endl; } else { i_output << "g" << i << "_i" << LEFT_ARRAY_SUBSCRIPT(output_type) << k + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << eq + 1 << ";" << endl; j_output << "g" << i << "_j" << LEFT_ARRAY_SUBSCRIPT(output_type) << k + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << col_idx + 1 << ";" << endl; v_output << "g" << i << "_v" << LEFT_ARRAY_SUBSCRIPT(output_type) << k + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << "="; d->writeOutput(v_output, output_type, temp_term_union, temporary_terms_idxs, tef_terms); v_output << ";" << endl; k++; } // Output symetric elements at order 2 if (i == 2 && vidx[1] != vidx[2]) { int col_idx_sym{getJacobianCol(vidx[2]) * getJacobianColsNbr() + getJacobianCol(vidx[1])}; if constexpr(isJuliaOutput(output_type)) d_output[2] << " g2[" << eq + 1 << "," << col_idx_sym + 1 << "] = " << "g2[" << eq + 1 << "," << col_idx + 1 << "]" << endl; else { i_output << "g" << i << "_i" << LEFT_ARRAY_SUBSCRIPT(output_type) << k + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << eq + 1 << ";" << endl; j_output << "g" << i << "_j" << LEFT_ARRAY_SUBSCRIPT(output_type) << k + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << col_idx_sym + 1 << ";" << endl; v_output << "g" << i << "_v" << LEFT_ARRAY_SUBSCRIPT(output_type) << k + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << "g" << i << "_v" << LEFT_ARRAY_SUBSCRIPT(output_type) << k-1 + ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << ";" << endl; k++; } } } if constexpr(!isJuliaOutput(output_type)) d_output[i] << i_output.str() << j_output.str() << v_output.str(); } if constexpr(isMatlabOutput(output_type)) { // Check that we don't have more than 32 nested parenthesis because MATLAB does not suppor this. See Issue #1201 map tmp_paren_vars; bool message_printed {false}; for (auto &it : tt_output) fixNestedParenthesis(it, tmp_paren_vars, message_printed); for (auto &it : d_output) fixNestedParenthesis(it, tmp_paren_vars, message_printed); } return { move(d_output), move(tt_output) }; } template void ModelTree::writePerBlockHelper(int blk, ostream &output, temporary_terms_t &temporary_terms) const { int block_recursive_size { blocks[blk].getRecursiveSize() }; // The equations deriv_node_temp_terms_t tef_terms; auto write_eq_tt = [&](int eq) { for (auto it : blocks_temporary_terms[blk][eq]) { if (dynamic_cast(it)) it->writeExternalFunctionOutput(output, output_type, temporary_terms, blocks_temporary_terms_idxs, tef_terms); output << " "; it->writeOutput(output, output_type, blocks_temporary_terms[blk][eq], blocks_temporary_terms_idxs, tef_terms); output << '='; it->writeOutput(output, output_type, temporary_terms, blocks_temporary_terms_idxs, tef_terms); temporary_terms.insert(it); output << ';' << endl; } }; for (int eq {0}; eq < blocks[blk].size; eq++) { write_eq_tt(eq); EquationType equ_type { getBlockEquationType(blk, eq) }; BinaryOpNode *e { getBlockEquationExpr(blk, eq) }; expr_t lhs { e->arg1 }, rhs { e->arg2 }; switch (blocks[blk].simulation_type) { case BlockSimulationType::evaluateBackward: case BlockSimulationType::evaluateForward: evaluation: if (equ_type == EquationType::evaluateRenormalized) { e = getBlockEquationRenormalizedExpr(blk, eq); lhs = e->arg1; rhs = e->arg2; } else if (equ_type != EquationType::evaluate) { cerr << "Type mismatch for equation " << getBlockEquationID(blk, eq)+1 << endl; exit(EXIT_FAILURE); } output << " "; lhs->writeOutput(output, output_type, temporary_terms, blocks_temporary_terms_idxs); output << '='; rhs->writeOutput(output, output_type, temporary_terms, blocks_temporary_terms_idxs); output << ';' << endl; break; case BlockSimulationType::solveBackwardSimple: case BlockSimulationType::solveForwardSimple: case BlockSimulationType::solveBackwardComplete: case BlockSimulationType::solveForwardComplete: case BlockSimulationType::solveTwoBoundariesComplete: case BlockSimulationType::solveTwoBoundariesSimple: if (eq < block_recursive_size) goto evaluation; output << " residual" << LEFT_ARRAY_SUBSCRIPT(output_type) << eq-block_recursive_size+ARRAY_SUBSCRIPT_OFFSET(output_type) << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=("; lhs->writeOutput(output, output_type, temporary_terms, blocks_temporary_terms_idxs); output << ")-("; rhs->writeOutput(output, output_type, temporary_terms, blocks_temporary_terms_idxs); output << ");" << endl; break; default: cerr << "Incorrect type for block " << blk+1 << endl; exit(EXIT_FAILURE); } } /* Write temporary terms for derivatives. This is done even for “evaluate” blocks, whose derivatives are not always computed at runtime; still those temporary terms may be needed by subsequent blocks. */ write_eq_tt(blocks[blk].size); } template tuple ModelTree::writeParamsDerivativesFileHelper() const { static_assert(!isCOutput(output_type), "C output is not implemented"); ostringstream tt_output; // Used for storing model temp vars and equations ostringstream rp_output; // 1st deriv. of residuals w.r.t. parameters ostringstream gp_output; // 1st deriv. of Jacobian w.r.t. parameters ostringstream rpp_output; // 2nd deriv of residuals w.r.t. parameters ostringstream gpp_output; // 2nd deriv of Jacobian w.r.t. parameters ostringstream hp_output; // 1st deriv. of Hessian w.r.t. parameters ostringstream g3p_output; // 1st deriv. of 3rd deriv. matrix w.r.t. parameters (only in dynamic case) temporary_terms_t temp_term_union; deriv_node_temp_terms_t tef_terms; writeModelLocalVariableTemporaryTerms(temp_term_union, params_derivs_temporary_terms_idxs, tt_output, tef_terms); for (const auto &[order, tts] : params_derivs_temporary_terms) writeTemporaryTerms(tts, temp_term_union, params_derivs_temporary_terms_idxs, tt_output, tef_terms); for (const auto &[indices, d1] : params_derivatives.find({ 0, 1 })->second) { auto [eq, param] { vectorToTuple<2>(indices) }; int param_col { getTypeSpecificIDByDerivID(param) + 1 }; rp_output << "rp" << LEFT_ARRAY_SUBSCRIPT(output_type) << eq+1 << ", " << param_col << RIGHT_ARRAY_SUBSCRIPT(output_type) << " = "; d1->writeOutput(rp_output, output_type, temp_term_union, params_derivs_temporary_terms_idxs, tef_terms); rp_output << ";" << endl; } for (const auto &[indices, d2] : params_derivatives.find({ 1, 1 })->second) { auto [eq, var, param] { vectorToTuple<3>(indices) }; int var_col { getJacobianCol(var) + 1 }; int param_col { getTypeSpecificIDByDerivID(param) + 1 }; gp_output << "gp" << LEFT_ARRAY_SUBSCRIPT(output_type) << eq+1 << ", " << var_col << ", " << param_col << RIGHT_ARRAY_SUBSCRIPT(output_type) << " = "; d2->writeOutput(gp_output, output_type, temp_term_union, params_derivs_temporary_terms_idxs, tef_terms); gp_output << ";" << endl; } for (int i {1}; const auto &[indices, d2] : params_derivatives.find({ 0, 2 })->second) { auto [eq, param1, param2] { vectorToTuple<3>(indices) }; int param1_col { getTypeSpecificIDByDerivID(param1) + 1 }; int param2_col { getTypeSpecificIDByDerivID(param2) + 1 }; rpp_output << "rpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",1" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << eq+1 << ";" << endl << "rpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",2" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param1_col << ";" << endl << "rpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",3" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param2_col << ";" << endl << "rpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",4" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "="; d2->writeOutput(rpp_output, output_type, temp_term_union, params_derivs_temporary_terms_idxs, tef_terms); rpp_output << ";" << endl; i++; if (param1 != param2) { // Treat symmetric elements rpp_output << "rpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",1" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << eq+1 << ";" << endl << "rpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",2" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param2_col << ";" << endl << "rpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",3" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param1_col << ";" << endl << "rpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",4" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=rpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i-1 << ",4" << RIGHT_ARRAY_SUBSCRIPT(output_type) << ";" << endl; i++; } } for (int i {1}; const auto &[indices, d2] : params_derivatives.find({ 1, 2 })->second) { auto [eq, var, param1, param2] { vectorToTuple<4>(indices) }; int var_col { getJacobianCol(var) + 1 }; int param1_col { getTypeSpecificIDByDerivID(param1) + 1 }; int param2_col { getTypeSpecificIDByDerivID(param2) + 1 }; gpp_output << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",1" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << eq+1 << ";" << endl << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",2" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << var_col << ";" << endl << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",3" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param1_col << ";" << endl << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",4" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param2_col << ";" << endl << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",5" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "="; d2->writeOutput(gpp_output, output_type, temp_term_union, params_derivs_temporary_terms_idxs, tef_terms); gpp_output << ";" << endl; i++; if (param1 != param2) { // Treat symmetric elements gpp_output << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",1" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << eq+1 << ";" << endl << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",2" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << var_col << ";" << endl << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",3" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param2_col << ";" << endl << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",4" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param1_col << ";" << endl << "gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",5" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=gpp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i-1 << ",5" << RIGHT_ARRAY_SUBSCRIPT(output_type) << ";" << endl; i++; } } for (int i {1}; const auto &[indices, d2] : params_derivatives.find({ 2, 1 })->second) { auto [eq, var1, var2, param] { vectorToTuple<4>(indices) }; int var1_col { getJacobianCol(var1) + 1 }; int var2_col { getJacobianCol(var2) + 1 }; int param_col { getTypeSpecificIDByDerivID(param) + 1 }; hp_output << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",1" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << eq+1 << ";" << endl << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",2" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << var1_col << ";" << endl << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",3" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << var2_col << ";" << endl << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",4" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param_col << ";" << endl << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",5" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "="; d2->writeOutput(hp_output, output_type, temp_term_union, params_derivs_temporary_terms_idxs, tef_terms); hp_output << ";" << endl; i++; if (var1 != var2) { // Treat symmetric elements hp_output << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",1" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << eq+1 << ";" << endl << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",2" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << var2_col << ";" << endl << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",3" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << var1_col << ";" << endl << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",4" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param_col << ";" << endl << "hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",5" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=hp" << LEFT_ARRAY_SUBSCRIPT(output_type) << i-1 << ",5" << RIGHT_ARRAY_SUBSCRIPT(output_type) << ";" << endl; i++; } } if constexpr(output_type == ExprNodeOutputType::matlabDynamicModel || output_type == ExprNodeOutputType::juliaDynamicModel) for (int i {1}; const auto &[indices, d2] : params_derivatives.find({ 3, 1 })->second) { auto [eq, var1, var2, var3, param] { vectorToTuple<5>(indices) }; int var1_col { getJacobianCol(var1) + 1 }; int var2_col { getJacobianCol(var2) + 1 }; int var3_col { getJacobianCol(var3) + 1 }; int param_col { getTypeSpecificIDByDerivID(param) + 1 }; g3p_output << "g3p" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",1" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << eq+1 << ";" << endl << "g3p" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",2" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << var1_col << ";" << endl << "g3p" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",3" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << var2_col << ";" << endl << "g3p" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",4" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << var3_col << ";" << endl << "g3p" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",5" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "=" << param_col << ";" << endl << "g3p" << LEFT_ARRAY_SUBSCRIPT(output_type) << i << ",6" << RIGHT_ARRAY_SUBSCRIPT(output_type) << "="; d2->writeOutput(g3p_output, output_type, temp_term_union, params_derivs_temporary_terms_idxs, tef_terms); g3p_output << ";" << endl; i++; } if constexpr(isMatlabOutput(output_type)) { // Check that we don't have more than 32 nested parenthesis because MATLAB does not support this. See Issue #1201 map tmp_paren_vars; bool message_printed {false}; fixNestedParenthesis(tt_output, tmp_paren_vars, message_printed); fixNestedParenthesis(rp_output, tmp_paren_vars, message_printed); fixNestedParenthesis(gp_output, tmp_paren_vars, message_printed); fixNestedParenthesis(rpp_output, tmp_paren_vars, message_printed); fixNestedParenthesis(gpp_output, tmp_paren_vars, message_printed); fixNestedParenthesis(hp_output, tmp_paren_vars, message_printed); fixNestedParenthesis(g3p_output, tmp_paren_vars, message_printed); } return { move(tt_output), move(rp_output), move(gp_output), move(rpp_output), move(gpp_output), move(hp_output), move(g3p_output) }; } template void ModelTree::writeBytecodeTemporaryTerms(const temporary_terms_t &tt, temporary_terms_t &temporary_terms_union, BytecodeWriter &code_file, deriv_node_temp_terms_t &tef_terms) const { for (auto it : tt) { if (dynamic_cast(it)) it->writeBytecodeExternalFunctionOutput(code_file, output_type, temporary_terms_union, temporary_terms_idxs, tef_terms); int idx {temporary_terms_idxs.at(it)}; code_file << FNUMEXPR_{ExpressionType::TemporaryTerm, idx}; it->writeBytecodeOutput(code_file, output_type, temporary_terms_union, temporary_terms_idxs, tef_terms); static_assert(output_type == ExprNodeBytecodeOutputType::dynamicModel || output_type == ExprNodeBytecodeOutputType::staticModel); if constexpr(output_type == ExprNodeBytecodeOutputType::dynamicModel) code_file << FSTPT_{idx}; else code_file << FSTPST_{idx}; temporary_terms_union.insert(it); } } template void ModelTree::writeBytecodeModelEquations(BytecodeWriter &code_file, const temporary_terms_t &temporary_terms, const deriv_node_temp_terms_t &tef_terms) const { for (int eq {0}; eq < static_cast(equations.size()); eq++) { BinaryOpNode *eq_node {equations[eq]}; expr_t lhs {eq_node->arg1}, rhs {eq_node->arg2}; code_file << FNUMEXPR_{ExpressionType::ModelEquation, eq}; // Test if the right hand side of the equation is empty. double vrhs {1.0}; try { vrhs = rhs->eval({}); } catch (ExprNode::EvalException &e) { } if (vrhs != 0) // The right hand side of the equation is not empty ⇒ residual=lhs-rhs { lhs->writeBytecodeOutput(code_file, output_type, temporary_terms, temporary_terms_idxs, tef_terms); rhs->writeBytecodeOutput(code_file, output_type, temporary_terms, temporary_terms_idxs, tef_terms); code_file << FBINARY_{BinaryOpcode::minus} << FSTPR_{eq}; } else // The right hand side of the equation is empty ⇒ residual=lhs { lhs->writeBytecodeOutput(code_file, output_type, temporary_terms, temporary_terms_idxs, tef_terms); code_file << FSTPR_{eq}; } } } template void ModelTree::writeBytecodeHelper(BytecodeWriter &code_file) const { constexpr ExprNodeBytecodeOutputType output_type { dynamic ? ExprNodeBytecodeOutputType::dynamicModel : ExprNodeBytecodeOutputType::staticModel }; temporary_terms_t temporary_terms_union; deriv_node_temp_terms_t tef_terms; writeBytecodeTemporaryTerms(temporary_terms_derivatives[0], temporary_terms_union, code_file, tef_terms); writeBytecodeModelEquations(code_file, temporary_terms_union, tef_terms); code_file << FENDEQU_{}; // Temporary terms for the Jacobian writeBytecodeTemporaryTerms(temporary_terms_derivatives[1], temporary_terms_union, code_file, tef_terms); // Get the current code_file position and jump if “evaluate” mode int pos_jmpifeval {code_file.getInstructionCounter()}; code_file << FJMPIFEVAL_{0}; // Use 0 as jump offset for the time being // The Jacobian in “simulate” mode vector>> my_derivatives(symbol_table.endo_nbr());; int count_u {symbol_table.endo_nbr()}; for (const auto &[indices, d1] : derivatives[1]) { auto [eq, deriv_id] {vectorToTuple<2>(indices)}; if (getTypeByDerivID(deriv_id) == SymbolType::endogenous) { int tsid {getTypeSpecificIDByDerivID(deriv_id)}; int lag {getLagByDerivID(deriv_id)}; if constexpr(dynamic) code_file << FNUMEXPR_{ExpressionType::FirstEndoDerivative, eq, tsid, lag}; else code_file << FNUMEXPR_{ExpressionType::FirstEndoDerivative, eq, tsid}; if (!my_derivatives[eq].size()) my_derivatives[eq].clear(); my_derivatives[eq].emplace_back(tsid, lag, count_u); d1->writeBytecodeOutput(code_file, output_type, temporary_terms_union, temporary_terms_idxs, tef_terms); if constexpr(dynamic) code_file << FSTPU_{count_u}; else code_file << FSTPSU_{count_u}; count_u++; } } for (int i {0}; i < symbol_table.endo_nbr(); i++) { code_file << FLDR_{i}; if (my_derivatives[i].size()) { for (bool first_term {true}; const auto &[tsid, lag, uidx] : my_derivatives[i]) { if constexpr(dynamic) code_file << FLDU_{uidx} << FLDV_{SymbolType::endogenous, tsid, lag}; else code_file << FLDSU_{uidx} << FLDSV_{SymbolType::endogenous, tsid}; code_file << FBINARY_{BinaryOpcode::times}; if (!exchange(first_term, false)) code_file << FBINARY_{BinaryOpcode::plus}; } code_file << FBINARY_{BinaryOpcode::minus}; } if constexpr(dynamic) code_file << FSTPU_{i}; else code_file << FSTPSU_{i}; } // Jump unconditionally after the block int pos_jmp {code_file.getInstructionCounter()}; code_file << FJMP_{0}; // Use 0 as jump offset for the time being // Update jump offset for previous JMPIFEVAL code_file.overwriteInstruction(pos_jmpifeval, FJMPIFEVAL_{pos_jmp-pos_jmpifeval}); // The Jacobian in “evaluate” mode for (const auto &[indices, d1] : derivatives[1]) { auto [eq, deriv_id] {vectorToTuple<2>(indices)}; int tsid {getTypeSpecificIDByDerivID(deriv_id)}; int lag {getLagByDerivID(deriv_id)}; SymbolType type {getTypeByDerivID(deriv_id)}; if constexpr(dynamic) { ExpressionType expr_type; switch (type) { case SymbolType::endogenous: expr_type = ExpressionType::FirstEndoDerivative; break; case SymbolType::exogenous: expr_type = ExpressionType::FirstExoDerivative; break; case SymbolType::exogenousDet: expr_type = ExpressionType::FirstExodetDerivative; break; default: assert(false); break; } code_file << FNUMEXPR_{expr_type, eq, tsid, lag}; } else { assert(type == SymbolType::endogenous); code_file << FNUMEXPR_{ExpressionType::FirstEndoDerivative, eq, tsid}; } d1->writeBytecodeOutput(code_file, output_type, temporary_terms_union, temporary_terms_idxs, tef_terms); if constexpr(dynamic) { // Bytecode MEX uses a separate matrix for exogenous and exodet Jacobians int jacob_col { type == SymbolType::endogenous ? getJacobianCol(deriv_id) : tsid }; code_file << FSTPG3_{eq, tsid, lag, jacob_col}; } else code_file << FSTPG2_{eq, tsid}; } // Update jump offset for previous JMP int pos_end_block {code_file.getInstructionCounter()}; code_file.overwriteInstruction(pos_jmp, FJMP_{pos_end_block-pos_jmp-1}); code_file << FENDBLOCK_{} << FEND_{}; } template void ModelTree::writeBlockBytecodeHelper(BytecodeWriter &code_file, int block) const { constexpr ExprNodeBytecodeOutputType output_type { dynamic ? ExprNodeBytecodeOutputType::dynamicModel : ExprNodeBytecodeOutputType::staticModel }; constexpr ExprNodeBytecodeOutputType assignment_lhs_output_type { dynamic ? ExprNodeBytecodeOutputType::dynamicAssignmentLHS : ExprNodeBytecodeOutputType::staticAssignmentLHS }; const BlockSimulationType simulation_type {blocks[block].simulation_type}; const int block_size {blocks[block].size}; const int block_mfs {blocks[block].mfs_size}; const int block_recursive {blocks[block].getRecursiveSize()}; temporary_terms_t temporary_terms_union; deriv_node_temp_terms_t tef_terms; auto write_eq_tt = [&](int eq) { for (auto it : blocks_temporary_terms[block][eq]) { if (dynamic_cast(it)) it->writeBytecodeExternalFunctionOutput(code_file, output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); code_file << FNUMEXPR_{ExpressionType::TemporaryTerm, blocks_temporary_terms_idxs.at(it)}; it->writeBytecodeOutput(code_file, output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); if constexpr(dynamic) code_file << FSTPT_{blocks_temporary_terms_idxs.at(it)}; else code_file << FSTPST_{blocks_temporary_terms_idxs.at(it)}; temporary_terms_union.insert(it); } }; // The equations for (int i {0}; i < block_size; i++) { write_eq_tt(i); switch (simulation_type) { evaluation: case BlockSimulationType::evaluateBackward: case BlockSimulationType::evaluateForward: code_file << FNUMEXPR_{ExpressionType::ModelEquation, getBlockEquationID(block, i)}; if (EquationType equ_type {getBlockEquationType(block, i)}; equ_type == EquationType::evaluate) { BinaryOpNode *eq_node {getBlockEquationExpr(block, i)}; expr_t lhs {eq_node->arg1}; expr_t rhs {eq_node->arg2}; rhs->writeBytecodeOutput(code_file, output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); lhs->writeBytecodeOutput(code_file, assignment_lhs_output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); } else if (equ_type == EquationType::evaluateRenormalized) { BinaryOpNode *eq_node {getBlockEquationRenormalizedExpr(block, i)}; expr_t lhs {eq_node->arg1}; expr_t rhs {eq_node->arg2}; rhs->writeBytecodeOutput(code_file, output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); lhs->writeBytecodeOutput(code_file, assignment_lhs_output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); } break; case BlockSimulationType::solveBackwardComplete: case BlockSimulationType::solveForwardComplete: case BlockSimulationType::solveTwoBoundariesComplete: case BlockSimulationType::solveTwoBoundariesSimple: if (i < block_recursive) goto evaluation; [[fallthrough]]; default: code_file << FNUMEXPR_{ExpressionType::ModelEquation, getBlockEquationID(block, i)}; BinaryOpNode *eq_node {getBlockEquationExpr(block, i)}; expr_t lhs {eq_node->arg1}; expr_t rhs {eq_node->arg2}; lhs->writeBytecodeOutput(code_file, output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); rhs->writeBytecodeOutput(code_file, output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); code_file << FBINARY_{BinaryOpcode::minus} << FSTPR_{i - block_recursive}; } } /* Write temporary terms for derivatives. This is done before FENDEQU, because residuals of a subsequent block may depend on temporary terms for the derivatives of the present block. Also note that in the case of “evaluate” blocks, derivatives are not computed in the “evaluate” mode; still their temporary terms must be computed even in that mode, because for the same reason as above they may be needed in subsequent blocks. */ write_eq_tt(blocks[block].size); code_file << FENDEQU_{}; // Get the current code_file position and jump if evaluating int pos_jmpifeval {code_file.getInstructionCounter()}; code_file << FJMPIFEVAL_{0}; // Use 0 as jump offset for the time being /* Write the derivatives for the “simulate” mode (not needed if the block is of type “evaluate backward/forward”) */ if (simulation_type != BlockSimulationType::evaluateBackward && simulation_type != BlockSimulationType::evaluateForward) { switch (simulation_type) { case BlockSimulationType::solveBackwardSimple: case BlockSimulationType::solveForwardSimple: { int eqr {getBlockEquationID(block, 0)}; int varr {getBlockVariableID(block, 0)}; code_file << FNUMEXPR_{ExpressionType::FirstEndoDerivative, eqr, varr, 0}; // Get contemporaneous derivative of the single variable in the block if (auto it { blocks_derivatives[block].find({ 0, 0, 0 }) }; it != blocks_derivatives[block].end()) it->second->writeBytecodeOutput(code_file, output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); else code_file << FLDZ_{}; code_file << FSTPG_{0}; } break; case BlockSimulationType::solveBackwardComplete: case BlockSimulationType::solveForwardComplete: case BlockSimulationType::solveTwoBoundariesComplete: case BlockSimulationType::solveTwoBoundariesSimple: { // For each equation, stores a list of tuples (index_u, var, lag) vector>> Uf(symbol_table.endo_nbr()); for (int count_u {block_mfs}; const auto &[indices, ignore] : blocks_derivatives[block]) { const auto &[eq, var, lag] {indices}; int eqr {getBlockEquationID(block, eq)}; int varr {getBlockVariableID(block, var)}; if (eq >= block_recursive && var >= block_recursive) { if constexpr(dynamic) if (lag != 0 && (simulation_type == BlockSimulationType::solveForwardComplete || simulation_type == BlockSimulationType::solveBackwardComplete)) continue; code_file << FNUMEXPR_{ExpressionType::FirstEndoDerivative, eqr, varr, lag}; if (auto it { blocks_derivatives[block].find({ eq, var, lag }) }; it != blocks_derivatives[block].end()) it->second->writeBytecodeOutput(code_file, output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); else code_file << FLDZ_{}; if constexpr(dynamic) code_file << FSTPU_{count_u}; else code_file << FSTPSU_{count_u}; Uf[eqr].emplace_back(count_u, varr, lag); count_u++; } } for (int i {0}; i < block_size; i++) if (i >= block_recursive) { code_file << FLDR_{i-block_recursive} << FLDZ_{}; int eqr {getBlockEquationID(block, i)}; for (const auto &[index_u, var, lag] : Uf[eqr]) { if constexpr(dynamic) code_file << FLDU_{index_u} << FLDV_{SymbolType::endogenous, var, lag}; else code_file << FLDSU_{index_u} << FLDSV_{SymbolType::endogenous, var}; code_file << FBINARY_{BinaryOpcode::times} << FBINARY_{BinaryOpcode::plus}; } code_file << FBINARY_{BinaryOpcode::minus}; if constexpr(dynamic) code_file << FSTPU_{i - block_recursive}; else code_file << FSTPSU_{i - block_recursive}; } } break; default: break; } } // Jump unconditionally after the block int pos_jmp {code_file.getInstructionCounter()}; code_file << FJMP_{0}; // Use 0 as jump offset for the time being // Update jump offset for previous JMPIFEVAL code_file.overwriteInstruction(pos_jmpifeval, FJMPIFEVAL_{pos_jmp-pos_jmpifeval}); // Write the derivatives for the “evaluate” mode for (const auto &[indices, d] : blocks_derivatives[block]) { const auto &[eq, var, lag] {indices}; int eqr {getBlockEquationID(block, eq)}; int varr {getBlockVariableID(block, var)}; code_file << FNUMEXPR_{ExpressionType::FirstEndoDerivative, eqr, varr, lag}; d->writeBytecodeOutput(code_file, output_type, temporary_terms_union, blocks_temporary_terms_idxs, tef_terms); if constexpr(dynamic) code_file << FSTPG3_{eq, var, lag, getBlockJacobianEndoCol(block, var, lag)}; else code_file << FSTPG2_{eq, getBlockJacobianEndoCol(block, var, lag)}; } /* Write derivatives w.r.t. exo, exo det and other endogenous, but only in dynamic mode */ writeBlockBytecodeAdditionalDerivatives(code_file, block, temporary_terms_union, tef_terms); // Update jump offset for previous JMP int pos_end_block {code_file.getInstructionCounter()}; code_file.overwriteInstruction(pos_jmp, FJMP_{pos_end_block-pos_jmp-1}); code_file << FENDBLOCK_{}; } template pair> ModelTree::writeJsonComputingPassOutputHelper(bool writeDetails) const { ostringstream mlv_output; // Used for storing model local vars vector d_output(derivatives.size()); // Derivatives output (at all orders, including 0=residual) deriv_node_temp_terms_t tef_terms; temporary_terms_t temp_term_union; writeJsonModelLocalVariables(mlv_output, true, tef_terms); writeJsonTemporaryTerms(temporary_terms_derivatives[0], temp_term_union, d_output[0], tef_terms, ""); d_output[0] << ", "; writeJsonModelEquations(d_output[0], true); int ncols { getJacobianColsNbr() }; for (size_t i {1}; i < derivatives.size(); i++) { string matrix_name { i == 1 ? "jacobian" : i == 2 ? "hessian" : i == 3 ? "third_derivative" : to_string(i) + "th_derivative"}; writeJsonTemporaryTerms(temporary_terms_derivatives[i], temp_term_union, d_output[i], tef_terms, matrix_name); temp_term_union.insert(temporary_terms_derivatives[i].begin(), temporary_terms_derivatives[i].end()); d_output[i] << R"(, ")" << matrix_name << R"(": {)" << R"( "nrows": )" << equations.size() << R"(, "ncols": )" << ncols << R"(, "entries": [)"; for (bool printed_something {false}; const auto &[vidx, d] : derivatives[i]) { if (exchange(printed_something, true)) d_output[i] << ", "; int eq { vidx[0] }; int col_idx {0}; for (size_t j {1}; j < vidx.size(); j++) { col_idx *= getJacobianColsNbr(); col_idx += getJacobianCol(vidx[j]); } if (writeDetails) d_output[i] << R"({"eq": )" << eq + 1; else d_output[i] << R"({"row": )" << eq + 1; d_output[i] << R"(, "col": )" << (i > 1 ? "[" : "") << col_idx + 1; if (i == 2 && vidx[1] != vidx[2]) // Symmetric elements in hessian { int col_idx_sym { getJacobianCol(vidx[2]) * getJacobianColsNbr() + getJacobianCol(vidx[1])}; d_output[i] << ", " << col_idx_sym + 1; } if (i > 1) d_output[i] << "]"; if (writeDetails) for (size_t j = 1; j < vidx.size(); j++) { d_output[i] << R"(, "var)" << (i > 1 ? to_string(j) : "") << R"(": ")" << getNameByDerivID(vidx[j]) << R"(")"; if constexpr(dynamic) d_output[i] << R"(, "shift)" << (i > 1 ? to_string(j) : "") << R"(": )" << getLagByDerivID(vidx[j]); } d_output[i] << R"(, "val": ")"; d->writeJsonOutput(d_output[i], temp_term_union, tef_terms); d_output[i] << R"("})" << endl; } d_output[i] << "]}"; ncols *= getJacobianColsNbr(); } return { move(mlv_output), move(d_output) }; } template tuple ModelTree::writeJsonParamsDerivativesHelper(bool writeDetails) const { ostringstream mlv_output; // Used for storing model local vars ostringstream tt_output; // Used for storing model temp vars and equations ostringstream rp_output; // 1st deriv. of residuals w.r.t. parameters ostringstream gp_output; // 1st deriv. of Jacobian w.r.t. parameters ostringstream rpp_output; // 2nd deriv of residuals w.r.t. parameters ostringstream gpp_output; // 2nd deriv of Jacobian w.r.t. parameters ostringstream hp_output; // 1st deriv. of Hessian w.r.t. parameters ostringstream g3p_output; // 1st deriv. of 3rd deriv. matrix w.r.t. parameters deriv_node_temp_terms_t tef_terms; writeJsonModelLocalVariables(mlv_output, true, tef_terms); temporary_terms_t temp_term_union; for (const auto &[order, tts] : params_derivs_temporary_terms) writeJsonTemporaryTerms(tts, temp_term_union, tt_output, tef_terms, "all"); rp_output << R"("deriv_wrt_params": {)" << R"( "neqs": )" << equations.size() << R"(, "nparamcols": )" << symbol_table.param_nbr() << R"(, "entries": [)"; for (bool printed_something {false}; const auto &[vidx, d] : params_derivatives.find({ 0, 1 })->second) { if (exchange(printed_something, true)) rp_output << ", "; auto [eq, param] { vectorToTuple<2>(vidx) }; int param_col { getTypeSpecificIDByDerivID(param) + 1 }; if (writeDetails) rp_output << R"({"eq": )" << eq + 1; else rp_output << R"({"row": )" << eq + 1; rp_output << R"(, "param_col": )" << param_col; if (writeDetails) rp_output << R"(, "param": ")" << getNameByDerivID(param) << R"(")"; rp_output << R"(, "val": ")"; d->writeJsonOutput(rp_output, temp_term_union, tef_terms); rp_output << R"("})" << endl; } rp_output << "]}"; gp_output << R"("deriv_jacobian_wrt_params": {)" << R"( "neqs": )" << equations.size() << R"(, "nvarcols": )" << getJacobianColsNbr() << R"(, "nparamcols": )" << symbol_table.param_nbr() << R"(, "entries": [)"; for (bool printed_something {false}; const auto &[vidx, d] : params_derivatives.find({ 1, 1 })->second) { if (exchange(printed_something, true)) gp_output << ", "; auto [eq, var, param] { vectorToTuple<3>(vidx) }; int var_col { getJacobianCol(var) + 1 }; int param_col { getTypeSpecificIDByDerivID(param) + 1 }; if (writeDetails) gp_output << R"({"eq": )" << eq + 1; else gp_output << R"({"row": )" << eq + 1; gp_output << R"(, "var_col": )" << var_col << R"(, "param_col": )" << param_col; if (writeDetails) { gp_output << R"(, "var": ")" << getNameByDerivID(var) << R"(")"; if constexpr(dynamic) gp_output << R"(, "lag": )" << getLagByDerivID(var); gp_output << R"(, "param": ")" << getNameByDerivID(param) << R"(")"; } gp_output << R"(, "val": ")"; d->writeJsonOutput(gp_output, temp_term_union, tef_terms); gp_output << R"("})" << endl; } gp_output << "]}"; rpp_output << R"("second_deriv_residuals_wrt_params": {)" << R"( "nrows": )" << equations.size() << R"(, "nparam1cols": )" << symbol_table.param_nbr() << R"(, "nparam2cols": )" << symbol_table.param_nbr() << R"(, "entries": [)"; for (bool printed_something {false}; const auto &[vidx, d] : params_derivatives.find({ 0, 2 })->second) { if (exchange(printed_something, true)) rpp_output << ", "; auto [eq, param1, param2] { vectorToTuple<3>(vidx) }; int param1_col { getTypeSpecificIDByDerivID(param1) + 1 }; int param2_col { getTypeSpecificIDByDerivID(param2) + 1 }; if (writeDetails) rpp_output << R"({"eq": )" << eq + 1; else rpp_output << R"({"row": )" << eq + 1; rpp_output << R"(, "param1_col": )" << param1_col << R"(, "param2_col": )" << param2_col; if (writeDetails) rpp_output << R"(, "param1": ")" << getNameByDerivID(param1) << R"(")" << R"(, "param2": ")" << getNameByDerivID(param2) << R"(")"; rpp_output << R"(, "val": ")"; d->writeJsonOutput(rpp_output, temp_term_union, tef_terms); rpp_output << R"("})" << endl; } rpp_output << "]}"; gpp_output << R"("second_deriv_jacobian_wrt_params": {)" << R"( "neqs": )" << equations.size() << R"(, "nvarcols": )" << getJacobianColsNbr() << R"(, "nparam1cols": )" << symbol_table.param_nbr() << R"(, "nparam2cols": )" << symbol_table.param_nbr() << R"(, "entries": [)"; for (bool printed_something {false}; const auto &[vidx, d] : params_derivatives.find({ 1, 2 })->second) { if (exchange(printed_something, true)) gpp_output << ", "; auto [eq, var, param1, param2] { vectorToTuple<4>(vidx) }; int var_col { getJacobianCol(var) + 1 }; int param1_col { getTypeSpecificIDByDerivID(param1) + 1 }; int param2_col { getTypeSpecificIDByDerivID(param2) + 1 }; if (writeDetails) gpp_output << R"({"eq": )" << eq + 1; else gpp_output << R"({"row": )" << eq + 1; gpp_output << R"(, "var_col": )" << var_col << R"(, "param1_col": )" << param1_col << R"(, "param2_col": )" << param2_col; if (writeDetails) { gpp_output << R"(, "var": ")" << getNameByDerivID(var) << R"(")"; if constexpr(dynamic) gpp_output << R"(, "lag": )" << getLagByDerivID(var); gpp_output << R"(, "param1": ")" << getNameByDerivID(param1) << R"(")" << R"(, "param2": ")" << getNameByDerivID(param2) << R"(")"; } gpp_output << R"(, "val": ")"; d->writeJsonOutput(gpp_output, temp_term_union, tef_terms); gpp_output << R"("})" << endl; } gpp_output << "]}" << endl; hp_output << R"("derivative_hessian_wrt_params": {)" << R"( "neqs": )" << equations.size() << R"(, "nvar1cols": )" << getJacobianColsNbr() << R"(, "nvar2cols": )" << getJacobianColsNbr() << R"(, "nparamcols": )" << symbol_table.param_nbr() << R"(, "entries": [)"; for (bool printed_something {false}; const auto &[vidx, d] : params_derivatives.find({ 2, 1 })->second) { if (exchange(printed_something, true)) hp_output << ", "; auto [eq, var1, var2, param] { vectorToTuple<4>(vidx) }; int var1_col { getJacobianCol(var1) + 1 }; int var2_col { getJacobianCol(var2) + 1 }; int param_col { getTypeSpecificIDByDerivID(param) + 1 }; if (writeDetails) hp_output << R"({"eq": )" << eq + 1; else hp_output << R"({"row": )" << eq + 1; hp_output << R"(, "var1_col": )" << var1_col << R"(, "var2_col": )" << var2_col << R"(, "param_col": )" << param_col; if (writeDetails) { hp_output << R"(, "var1": ")" << getNameByDerivID(var1) << R"(")"; if constexpr(dynamic) hp_output << R"(, "lag1": )" << getLagByDerivID(var1); hp_output << R"(, "var2": ")" << getNameByDerivID(var2) << R"(")"; if constexpr(dynamic) hp_output << R"(, "lag2": )" << getLagByDerivID(var2); hp_output << R"(, "param": ")" << getNameByDerivID(param) << R"(")"; } hp_output << R"(, "val": ")"; d->writeJsonOutput(hp_output, temp_term_union, tef_terms); hp_output << R"("})" << endl; } hp_output << "]}" << endl; if constexpr(dynamic) { g3p_output << R"("derivative_g3_wrt_params": {)" << R"( "neqs": )" << equations.size() << R"(, "nvar1cols": )" << getJacobianColsNbr() << R"(, "nvar2cols": )" << getJacobianColsNbr() << R"(, "nvar3cols": )" << getJacobianColsNbr() << R"(, "nparamcols": )" << symbol_table.param_nbr() << R"(, "entries": [)"; for (bool printed_something {false}; const auto &[vidx, d] : params_derivatives.find({ 3, 1 })->second) { if (exchange(printed_something, true)) g3p_output << ", "; auto [eq, var1, var2, var3, param] { vectorToTuple<5>(vidx) }; int var1_col { getJacobianCol(var1) + 1 }; int var2_col { getJacobianCol(var2) + 1 }; int var3_col { getJacobianCol(var3) + 1 }; int param_col { getTypeSpecificIDByDerivID(param) + 1 }; if (writeDetails) g3p_output << R"({"eq": )" << eq + 1; else g3p_output << R"({"row": )" << eq + 1; g3p_output << R"(, "var1_col": )" << var1_col + 1 << R"(, "var2_col": )" << var2_col + 1 << R"(, "var3_col": )" << var3_col + 1 << R"(, "param_col": )" << param_col + 1; if (writeDetails) g3p_output << R"(, "var1": ")" << getNameByDerivID(var1) << R"(")" << R"(, "lag1": )" << getLagByDerivID(var1) << R"(, "var2": ")" << getNameByDerivID(var2) << R"(")" << R"(, "lag2": )" << getLagByDerivID(var2) << R"(, "var3": ")" << getNameByDerivID(var3) << R"(")" << R"(, "lag3": )" << getLagByDerivID(var3) << R"(, "param": ")" << getNameByDerivID(param) << R"(")"; g3p_output << R"(, "val": ")"; d->writeJsonOutput(g3p_output, temp_term_union, tef_terms); g3p_output << R"("})" << endl; } g3p_output << "]}" << endl; } return { move(mlv_output), move(tt_output), move(rp_output), move(gp_output), move(rpp_output), move(gpp_output), move(hp_output), move(g3p_output) }; } #endif