model_c_source_gen.hpp

#ifndef CPPAD_CG_MODEL_C_SOURCE_GEN_INCLUDED
#define CPPAD_CG_MODEL_C_SOURCE_GEN_INCLUDED
/* --------------------------------------------------------------------------
 *  CppADCodeGen: C++ Algorithmic Differentiation with Source Code Generation:
 *    Copyright (C) 2012 Ciengis
 *    Copyright (C) 2019 Joao Leal
 *
 *  CppADCodeGen is distributed under multiple licenses:
 *
 *   - Eclipse Public License Version 1.0 (EPL1), and
 *   - GNU General Public License Version 3 (GPL3).
 *
 *  EPL1 terms and conditions can be found in the file "epl-v10.txt", while
 *  terms and conditions for the GPL3 can be found in the file "gpl3.txt".
 * ----------------------------------------------------------------------------
 * Author: Joao Leal
 */

namespace CppAD {
namespace cg {

class CompressedVectorInfo {
public:
    std::vector<std::set<size_t> > locations; // location of each index in the compressed array
    std::vector<size_t> indexes; // row or column index
    bool ordered;
};

template<class Base>
class ModelCSourceGen {
    using CGBase = CppAD::cg::CG<Base>;
    using ADCG = CppAD::AD<CGBase>;
    using SparsitySetType = std::vector<std::set<size_t> >;
    using TapeVarType = std::pair<size_t, size_t>; // tape independent -> reference orig independent (temporaries only)
public:
    static const std::string FUNCTION_FORWAD_ZERO;
    static const std::string FUNCTION_JACOBIAN;
    static const std::string FUNCTION_HESSIAN;
    static const std::string FUNCTION_FORWARD_ONE;
    static const std::string FUNCTION_REVERSE_ONE;
    static const std::string FUNCTION_REVERSE_TWO;
    static const std::string FUNCTION_SPARSE_JACOBIAN;
    static const std::string FUNCTION_SPARSE_HESSIAN;
    static const std::string FUNCTION_JACOBIAN_SPARSITY;
    static const std::string FUNCTION_HESSIAN_SPARSITY;
    static const std::string FUNCTION_HESSIAN_SPARSITY2;
    static const std::string FUNCTION_SPARSE_FORWARD_ONE;
    static const std::string FUNCTION_SPARSE_REVERSE_ONE;
    static const std::string FUNCTION_SPARSE_REVERSE_TWO;
    static const std::string FUNCTION_FORWARD_ONE_SPARSITY;
    static const std::string FUNCTION_REVERSE_ONE_SPARSITY;
    static const std::string FUNCTION_REVERSE_TWO_SPARSITY;
    static const std::string FUNCTION_INFO;
    static const std::string FUNCTION_ATOMIC_FUNC_NAMES;
protected:
    static const std::string CONST;

    class Position {
    public:
        bool defined;
        std::vector<size_t> row;
        std::vector<size_t> col;

        inline Position() :
            defined(false) {
        }

        inline Position(const std::vector<size_t>& r, const std::vector<size_t>& c) :
            defined(true),
            row(r),
            col(c) {
            CPPADCG_ASSERT_KNOWN(r.size() == c.size(), "The number of row indexes must be the same as the number of column indexes.")
        }

        template<class VectorSet>
        inline Position(const VectorSet& elements) :
            defined(true) {
            size_t nnz = 0;
            for (size_t i = 0; i < elements.size(); i++) {
                nnz += elements[i].size();
            }
            row.resize(nnz);
            col.resize(nnz);

            nnz = 0;
            for (size_t i = 0; i < elements.size(); i++) {
                for (size_t it : elements[i]) {
                    row[nnz] = i;
                    col[nnz] = it;
                    nnz++;
                }
            }
        }
    };

    class LocalSparsityInfo {
    public:
        SparsitySetType sparsity;
        // rows (in a custom order)
        std::vector<size_t> rows;
        // columns (in a custom order)
        std::vector<size_t> cols;
    };

    class Color {
    public:
        std::set<size_t> rows;
        std::map<size_t, size_t> column2Row;
        std::map<size_t, std::set<size_t> > row2Columns;
        std::set<size_t> forbiddenRows;
    };

protected:
    ADFun<CGBase>& _fun;
    LoopFreeModel<Base>* _funNoLoops;
    std::set<LoopModel<Base>*> _loopTapes;
    std::string _name;
    const std::string _baseTypeName;
    size_t _parameterPrecision;
    std::vector<Base> _x;
    bool _multiThreading;
    bool _zero;
    bool _zeroEvaluated;
    bool _jacobian;
    bool _hessian;
    bool _sparseJacobian;
    bool _sparseHessian;
    bool _hessianByEquation;
    bool _forwardOne;
    bool _reverseOne;
    bool _reverseTwo;
    bool _sparseJacobianReusesOne;
    bool _sparseHessianReusesRev2;
    JacobianADMode _jacMode;
    Position _custom_jac;
    LocalSparsityInfo _jacSparsity;
    Position _custom_hess;
    LocalSparsityInfo _hessSparsity;
    std::vector<LocalSparsityInfo> _hessSparsities;
    std::vector<std::string> _atomicFunctions;
    std::map<size_t, AtomicUseInfo<Base> >* _atomicsInfo;
    std::ostringstream _cache;
    size_t _maxAssignPerFunc;
    size_t _maxOperationsPerAssignment;
    std::vector<std::set<size_t> > _relatedDepCandidates;
    std::map<LoopModel<Base>*, std::map<size_t, std::map<size_t, std::set<size_t> > > > _loopFor1Groups;
    std::map<size_t, std::set<size_t> > _nonLoopFor1Elements;
    std::map<LoopModel<Base>*, std::map<size_t, std::map<size_t, std::set<size_t> > > > _loopRev1Groups;
    std::map<size_t, std::set<size_t> > _nonLoopRev1Elements;
    std::map<LoopModel<Base>*, std::map<size_t, std::map<size_t, std::set<size_t> > > > _loopRev2Groups;
    std::map<size_t, std::set<size_t> > _nonLoopRev2Elements;
    JobTimer* _jobTimer;
    std::map<std::string, std::string> _sources;
public:

    ModelCSourceGen(ADFun<CppAD::cg::CG<Base> >& fun,
                    std::string model) :
        _fun(fun),
        _funNoLoops(nullptr),
        _name(std::move(model)),
        _baseTypeName(ModelCSourceGen<Base>::baseTypeName()),
        _parameterPrecision(std::numeric_limits<Base>::digits10),
        _multiThreading(true),
        _zero(true),
        _zeroEvaluated(false),
        _jacobian(false),
        _hessian(false),
        _sparseJacobian(false),
        _sparseHessian(false),
        _hessianByEquation(false),
        _forwardOne(false),
        _reverseOne(false),
        _reverseTwo(false),
        _sparseJacobianReusesOne(true),
        _sparseHessianReusesRev2(true),
        _jacMode(JacobianADMode::Automatic),
        _atomicsInfo(nullptr),
        _maxAssignPerFunc(20000),
        _maxOperationsPerAssignment(1000),
        _jobTimer(nullptr) {

        CPPADCG_ASSERT_KNOWN(!_name.empty(), "Model name cannot be empty")
        CPPADCG_ASSERT_KNOWN((_name[0] >= 'a' && _name[0] <= 'z') ||
                             (_name[0] >= 'A' && _name[0] <= 'Z'),
                             "Invalid model name character")
        for (size_t i = 1; i < _name.size(); i++) {
            char c = _name[i];
            CPPADCG_ASSERT_KNOWN((c >= 'a' && c <= 'z') ||
                                 (c >= 'A' && c <= 'Z') ||
                                 (c >= '0' && c <= '9') ||
                                 c == '_'
                                 , "Invalid model name character")
        }
    }

    ModelCSourceGen(const ModelCSourceGen&) = delete;
    ModelCSourceGen& operator=(const ModelCSourceGen&) = delete;

    inline const std::string& getName() const {
        return _name;
    }

    template<class VectorBase>
    inline void setTypicalIndependentValues(const VectorBase& x) {
        CPPAD_ASSERT_KNOWN(x.size() == 0 || x.size() == _fun.Domain(),
                           "Invalid independent variable vector size")
        _x.resize(x.size());
        for (size_t i = 0; i < x.size(); i++) {
            _x[i] = x[i];
        }
    }

    inline void setRelatedDependents(const std::vector<std::set<size_t> >& relatedDepCandidates) {
        _relatedDepCandidates = relatedDepCandidates;
    }

    inline const std::vector<std::set<size_t> >& getRelatedDependents() const {
        return _relatedDepCandidates;
    }

    virtual size_t getParameterPrecision() const {
        return _parameterPrecision;
    }

    virtual void setParameterPrecision(size_t p) {
        _parameterPrecision = p;
    }

    inline bool isMultiThreading() const {
        return _multiThreading;
    }

    inline void setMultiThreading(bool multiThreading) {
        _multiThreading = multiThreading;
    }

    inline bool isJacobianMultiThreadingEnabled() const {
        return _multiThreading && _loopTapes.empty() && _sparseJacobian && _sparseJacobianReusesOne && (_forwardOne || _reverseOne);
    }

    inline bool isHessianMultiThreadingEnabled() const {
        return _multiThreading && _loopTapes.empty() && _sparseHessian && _sparseHessianReusesRev2 && _reverseTwo;
    }

    inline bool isCreateHessian() const {
        return _hessian;
    }

    inline void setCreateHessian(bool create) {
        _hessian = create;
    }

    inline JacobianADMode getJacobianADMode() const {
        return _jacMode;
    }

    inline void setJacobianADMode(JacobianADMode mode) {
        _jacMode = mode;
    }

    inline bool isCreateJacobian() const {
        return _jacobian;
    }

    inline void setCreateJacobian(bool create) {
        _jacobian = create;
    }

    inline bool isCreateSparseHessian() const {
        return _sparseHessian;
    }

    inline void setCreateSparseHessian(bool create) {
        _sparseHessian = create;
    }

    inline bool isSparseHessianReusesRev2() const {
        return _sparseHessianReusesRev2;
    }

    inline void setSparseHessianReusesRev2(bool reuse) {
        _sparseHessianReusesRev2 = reuse;
    }

    inline bool isCreateHessianSparsityByEquation() const {
        return _hessianByEquation;
    }

    inline void setCreateHessianSparsityByEquation(bool create) {
        _hessianByEquation = create;
    }

    inline bool isCreateSparseJacobian() const {
        return _sparseJacobian;
    }

    inline void setCreateSparseJacobian(bool create) {
        _sparseJacobian = create;
    }

    inline bool isSparseJacobianReuse1stOrderPasses() const {
        return _sparseJacobianReusesOne;
    }

    inline void setSparseJacobianReuse1stOrderPasses(bool reuse) {
        _sparseJacobianReusesOne = reuse;
    }

    inline bool isCreateForwardZero() const {
        return _zero;
    }

    inline void setCreateForwardZero(bool create) {
        _zero = create;
    }

    inline bool isCreateSparseForwardOne() const {
        return _forwardOne;
    }

    inline void setCreateForwardOne(bool create) {
        _forwardOne = create;
    }

    inline bool isCreateReverseOne() const {
        return _reverseOne;
    }

    inline void setCreateReverseOne(bool create) {
        _reverseOne = create;
    }

    inline bool isCreateReverseTwo() const {
        return _reverseOne;
    }

    inline void setCreateReverseTwo(bool create) {
        _reverseTwo = create;
    }

    inline void setCustomSparseJacobianElements(const std::vector<size_t>& row,
                                                const std::vector<size_t>& col) {
        _custom_jac = Position(row, col);
    }

    template<class VectorSet>
    inline void setCustomSparseJacobianElements(const VectorSet& elements) {
        _custom_jac = Position(elements);
    }

    inline void setCustomSparseHessianElements(const std::vector<size_t>& row,
                                               const std::vector<size_t>& col) {
        _custom_hess = Position(row, col);
    }

    template<class VectorSet>
    inline void setCustomSparseHessianElements(const VectorSet& elements) {
        _custom_hess = Position(elements);
    }

    inline size_t getMaxAssignmentsPerFunc() const {
        return _maxAssignPerFunc;
    }

    inline void setMaxAssignmentsPerFunc(size_t maxAssignPerFunc) {
        _maxAssignPerFunc = maxAssignPerFunc;
    }

    inline size_t getMaxOperationsPerAssignment() const {
        return _maxOperationsPerAssignment;
    }

    inline void setMaxOperationsPerAssignment(size_t maxOperationsPerAssignment) {
        _maxOperationsPerAssignment = maxOperationsPerAssignment;
    }

    inline virtual ~ModelCSourceGen() {
        delete _funNoLoops;
        delete _atomicsInfo;

        for (LoopModel<Base>* it : _loopTapes) {
            delete it;
        }
    }

public:
    static inline std::string baseTypeName();

    template<class T>
    static void generateFunctionDeclarationSource(std::ostringstream& cache,
                                                  const std::string& model_function,
                                                  const std::string& suffix,
                                                  const std::map<size_t, T>& elements,
                                                  const std::string& argsDcl);

protected:

    virtual VariableNameGenerator<Base>* createVariableNameGenerator(const std::string& depName = "y",
                                                                     const std::string& indepName = "x",
                                                                     const std::string& tmpName = "v",
                                                                     const std::string& tmpArrayName = "array");

    const std::map<std::string, std::string>& getSources(MultiThreadingType multiThreadingType,
                                                         JobTimer* timer);

    virtual void generateSources(MultiThreadingType multiThreadingType,
                                 JobTimer* timer = nullptr);

    virtual void generateLoops();

    virtual void generateInfoSource();

    virtual void generateAtomicFuncNames();

    virtual bool isAtomicsUsed();

    virtual const std::map<size_t, AtomicUseInfo<Base> >& getAtomicsInfo();

    /***********************************************************************
     * zero order (the original model)
     **********************************************************************/

    virtual void generateZeroSource();

    virtual std::vector<CGBase> prepareForward0WithLoops(CodeHandler<Base>& handler,
                                                         const std::vector<CGBase>& x);

    /***********************************************************************
     * Jacobian
     **********************************************************************/

    virtual void generateJacobianSource();

    virtual void generateSparseJacobianSource(MultiThreadingType multiThreadingType);

    virtual void generateSparseJacobianSource(bool forward);

    virtual void generateSparseJacobianForRevSource(bool forward,
                                                    MultiThreadingType multiThreadingType);

    virtual std::string generateSparseJacobianForRevSingleThreadSource(const std::string& functionName,
                                                                       std::map<size_t, CompressedVectorInfo> jacInfo,
                                                                       size_t maxCompressedSize,
                                                                       const std::string& functionRevFor,
                                                                       const std::string& revForSuffix,
                                                                       bool forward);

    virtual std::string generateSparseJacobianForRevMultiThreadSource(const std::string& functionName,
                                                                      std::map<size_t, CompressedVectorInfo> jacInfo,
                                                                      size_t maxCompressedSize,
                                                                      const std::string& functionRevFor,
                                                                      const std::string& revForSuffix,
                                                                      bool forward,
                                                                      MultiThreadingType multiThreadingType);
    virtual std::vector<CGBase> prepareSparseJacobianWithLoops(CodeHandler<Base>& handler,
                                                               const std::vector<CGBase>& x,
                                                               bool forward);

    inline void prepareSparseJacobianRowWithLoops(CodeHandler<Base>& handler,
                                                  LoopModel<Base>& lModel,
                                                  size_t tapeI,
                                                  const loops::JacobianWithLoopsRowInfo& rowInfo,
                                                  const std::vector<std::map<size_t, CGBase> >& dyiDxtape,
                                                  const std::vector<std::map<size_t, CGBase> >& dzDx,
                                                  const CGBase& py,
                                                  IndexOperationNode<Base>& iterationIndexOp,
                                                  std::vector<loops::IfElseInfo<Base> >& ifElses,
                                                  size_t& jacLE,
                                                  std::vector<std::pair<CG<Base>, IndexPattern*> >& indexedLoopResults,
                                                  std::set<size_t>& allLocations);

    inline void analyseSparseJacobianWithLoops(const std::vector<size_t>& rows,
                                               const std::vector<size_t>& cols,
                                               const std::vector<size_t>& location,
                                               SparsitySetType& noLoopEvalSparsity,
                                               std::vector<std::map<size_t, std::set<size_t> > >& noLoopEvalLocations,
                                               std::map<LoopModel<Base>*, SparsitySetType>& loopsEvalSparsities,
                                               std::map<LoopModel<Base>*, std::vector<loops::JacobianWithLoopsRowInfo> >& loopEqInfo);


    virtual void generateSparseJacobianWithLoopsSourceFromForRev(const std::map<size_t, CompressedVectorInfo>& jacInfo,
                                                                 size_t maxCompressedSize,
                                                                 const std::string& localFunctionTypeName,
                                                                 const std::string& suffix,
                                                                 const std::string& keyName,
                                                                 const std::map<size_t, std::set<size_t> >& nonLoopElements,
                                                                 const std::map<LoopModel<Base>*, std::map<size_t, std::map<size_t, std::set<size_t> > > >& loopGroups,
                                                                 void (*generateLocalFunctionName)(std::ostringstream& cache, const std::string& modelName, const LoopModel<Base>& loop, size_t g));

    inline virtual void generateFunctionNameLoopFor1(std::ostringstream& cache,
                                                     const LoopModel<Base>& loop,
                                                     size_t g);

    inline static void generateFunctionNameLoopFor1(std::ostringstream& cache,
                                                    const std::string& modelName,
                                                    const LoopModel<Base>& loop,
                                                    size_t g);

    inline virtual void generateFunctionNameLoopRev1(std::ostringstream& cache,
                                                     const LoopModel<Base>& loop,
                                                     size_t i);

    inline static void generateFunctionNameLoopRev1(std::ostringstream& cache,
                                                    const std::string& modelName,
                                                    const LoopModel<Base>& loop,
                                                    size_t i);

    /***********************************************************************
     * Hessian
     **********************************************************************/

    virtual void generateHessianSource();

    virtual void generateSparseHessianSource(MultiThreadingType multiThreadingType);

    virtual void generateSparseHessianSourceDirectly();

    virtual void generateSparseHessianSourceFromRev2(MultiThreadingType multiThreadingType);

    virtual std::string generateSparseHessianRev2SingleThreadSource(const std::string& functionName,
                                                                    std::map<size_t, CompressedVectorInfo> hessInfo,
                                                                    size_t maxCompressedSize,
                                                                    const std::string& functionRev2,
                                                                    const std::string& rev2Suffix);

    virtual std::string generateSparseHessianRev2MultiThreadSource(const std::string& functionName,
                                                                   std::map<size_t, CompressedVectorInfo> hessInfo,
                                                                   size_t maxCompressedSize,
                                                                   const std::string& functionRev2,
                                                                   const std::string& rev2Suffix,
                                                                   MultiThreadingType multiThreadingType);

    virtual void determineSecondOrderElements4Eval(std::vector<size_t>& userRows,
                                                   std::vector<size_t>& userCols);

    virtual std::vector<CGBase> prepareSparseHessianWithLoops(CodeHandler<Base>& handler,
                                                              std::vector<CGBase>& indVars,
                                                              std::vector<CGBase>& w,
                                                              const std::vector<size_t>& lowerHessRows,
                                                              const std::vector<size_t>& lowerHessCols,
                                                              const std::vector<size_t>& lowerHessOrder,
                                                              const std::map<size_t, size_t>& duplicates);

    inline void analyseSparseHessianWithLoops(const std::vector<size_t>& lowerHessRows,
                                              const std::vector<size_t>& lowerHessCols,
                                              const std::vector<size_t>& lowerHessOrder,
                                              SparsitySetType& noLoopEvalJacSparsity,
                                              SparsitySetType& noLoopEvalHessSparsity,
                                              std::vector<std::map<size_t, std::set<size_t> > >& noLoopEvalHessLocations,
                                              std::map<LoopModel<Base>*, loops::HessianWithLoopsInfo<Base> >& loopHessInfo,
                                              bool useSymmetry);

    inline virtual void generateSparseHessianWithLoopsSourceFromRev2(const std::map<size_t, CompressedVectorInfo>& hessInfo,
                                                                     size_t maxCompressedSize);

    inline virtual void generateFunctionNameLoopRev2(std::ostringstream& cache,
                                                     const LoopModel<Base>& loop,
                                                     size_t g);

    static inline void generateFunctionNameLoopRev2(std::ostringstream& cache,
                                                    const std::string& modelName,
                                                    const LoopModel<Base>& loop,
                                                    size_t g);

    /***********************************************************************
     * Sparsities for forward/reverse
     **********************************************************************/

    virtual void generateSparsity1DSource(const std::string& function,
                                          const std::vector<size_t>& sparsity);

    virtual void generateSparsity2DSource(const std::string& function,
                                          const LocalSparsityInfo& sparsity);

    virtual void generateSparsity2DSource2(const std::string& function,
                                           const std::vector<LocalSparsityInfo>& sparsities);

    virtual void generateSparsity1DSource2(const std::string& function,
                                           const std::map<size_t, std::vector<size_t> >& rows);

    /***********************************************************************
     * Forward 1 mode
     **********************************************************************/

    virtual void generateSparseForwardOneSources();

    virtual void generateSparseForwardOneSourcesWithAtomics(const std::map<size_t, std::vector<size_t> >& elements);

    virtual void generateSparseForwardOneSourcesNoAtomics(const std::map<size_t, std::vector<size_t> >& elements);

    virtual void generateForwardOneSources();

    virtual void prepareSparseForwardOneWithLoops(const std::map<size_t, std::vector<size_t> >& elements);

    virtual void createForwardOneWithLoopsNL(CodeHandler<Base>& handler,
                                             size_t j,
                                             std::vector<CG<Base> >& jacCol);


    inline static std::map<size_t, std::map<size_t, CG<Base> > > generateLoopFor1Jac(ADFun<CGBase>& fun,
                                                                                     const SparsitySetType& sparsity,
                                                                                     const SparsitySetType& evalSparsity,
                                                                                     const std::vector<CGBase>& xl,
                                                                                     bool constainsAtomics);

    /***********************************************************************
     * Reverse 1 mode
     **********************************************************************/

    virtual void generateSparseReverseOneSources();

    virtual void generateSparseReverseOneSourcesWithAtomics(const std::map<size_t, std::vector<size_t> >& elements);

    virtual void generateSparseReverseOneSourcesNoAtomics(const std::map<size_t, std::vector<size_t> >& elements);

    virtual void generateReverseOneSources();

    virtual void prepareSparseReverseOneWithLoops(const std::map<size_t, std::vector<size_t> >& elements);

    virtual void createReverseOneWithLoopsNL(CodeHandler<Base>& handler,
                                             size_t i,
                                             std::vector<CG<Base> >& jacRow);

    inline static std::vector<std::map<size_t, CGBase> > generateLoopRev1Jac(ADFun<CGBase>& fun,
                                                                             const SparsitySetType& sparsity,
                                                                             const SparsitySetType& evalSparsity,
                                                                             const std::vector<CGBase>& xl,
                                                                             bool constainsAtomics);

    /***********************************************************************
     * Reverse 2 mode
     **********************************************************************/

    virtual void generateSparseReverseTwoSources();

    virtual void generateSparseReverseTwoSourcesWithAtomics(const std::map<size_t, std::vector<size_t> >& elements);

    virtual void generateSparseReverseTwoSourcesNoAtomics(const std::map<size_t, std::vector<size_t> >& elements,
                                                          const std::vector<size_t>& evalRows,
                                                          const std::vector<size_t>& evalCols);

    virtual void generateReverseTwoSources();

    virtual void generateGlobalDirectionalFunctionSource(const std::string& function,
                                                         const std::string& function2_suffix,
                                                         const std::string& function_sparsity,
                                                         const std::map<size_t, std::vector<size_t> >& elements);

    virtual void prepareSparseReverseTwoWithLoops(const std::map<size_t, std::vector<size_t> >& elements);

    /***********************************************************************
     * Sparsities
     **********************************************************************/

    virtual void determineJacobianSparsity();

    virtual void generateJacobianSparsitySource();

    virtual void determineHessianSparsity();

    inline std::vector<ModelCSourceGen<Base>::Color> colorByRow(const std::set<size_t>& columns,
                                                                const SparsitySetType& sparsity);

    virtual void generateHessianSparsitySource();


    static inline std::map<size_t, std::vector<std::set<size_t> > > determineOrderByCol(const std::map<size_t, std::vector<size_t> >& elements,
                                                                                        const LocalSparsityInfo& sparsity);

    static inline std::map<size_t, std::vector<std::set<size_t> > > determineOrderByCol(const std::map<size_t, std::vector<size_t> >& elements,
                                                                                        const std::vector<size_t>& userRows,
                                                                                        const std::vector<size_t>& userCols);

    static inline std::vector<std::set<size_t> > determineOrderByCol(size_t col,
                                                                     const std::vector<size_t>& colElements,
                                                                     const std::vector<size_t>& userRows,
                                                                     const std::vector<size_t>& userCols);

    static inline std::map<size_t, std::vector<std::set<size_t> > > determineOrderByRow(const std::map<size_t, std::vector<size_t> >& elements,
                                                                                        const LocalSparsityInfo& sparsity);

    static inline std::map<size_t, std::vector<std::set<size_t> > > determineOrderByRow(const std::map<size_t, std::vector<size_t> >& elements,
                                                                                        const std::vector<size_t>& userRows,
                                                                                        const std::vector<size_t>& userCols);

    static inline std::vector<std::set<size_t> > determineOrderByRow(size_t row,
                                                                     const std::vector<size_t>& rowsElements,
                                                                     const std::vector<size_t>& userRows,
                                                                     const std::vector<size_t>& userCols);

    /***********************************************************************
     * Multi-threading
     **********************************************************************/

    static void printFileStartPThreads(std::ostringstream& cache,
                                       const std::string& baseTypeName);

    static void printFunctionStartPThreads(std::ostringstream& cache,
                                           size_t size);

    static void printFunctionEndPThreads(std::ostringstream& cache,
                                         size_t size);

    static void printFileStartOpenMP(std::ostringstream& cache);

    static void printFunctionStartOpenMP(std::ostringstream& cache,
                                         size_t size);

    static void printLoopStartOpenMP(std::ostringstream& cache,
                                     size_t size);

    static void printLoopEndOpenMP(std::ostringstream& cache,
                                   size_t size);

    inline void startingJob(const std::string& jobName,
                            const JobType& type = JobTypeHolder<>::DEFAULT);

    inline void finishedJob();

    friend class
    ModelLibraryCSourceGen<Base>;

    friend class
    ModelLibraryProcessor<Base>;
};

} // END cg namespace
} // END CppAD namespace

#endif