function_internal.hpp
1 /*
2  * This file is part of CasADi.
3  *
4  * CasADi -- A symbolic framework for dynamic optimization.
5  * Copyright (C) 2010-2014 Joel Andersson, Joris Gillis, Moritz Diehl,
6  * KU Leuven. All rights reserved.
7  * Copyright (C) 2011-2014 Greg Horn
8  *
9  * CasADi is free software; you can redistribute it and/or
10  * modify it under the terms of the GNU Lesser General Public
11  * License as published by the Free Software Foundation; either
12  * version 3 of the License, or (at your option) any later version.
13  *
14  * CasADi is distributed in the hope that it will be useful,
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17  * Lesser General Public License for more details.
18  *
19  * You should have received a copy of the GNU Lesser General Public
20  * License along with CasADi; if not, write to the Free Software
21  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
22  *
23  */
24 
25 
26 #ifndef CASADI_FUNCTION_INTERNAL_HPP
27 #define CASADI_FUNCTION_INTERNAL_HPP
28 
29 #include "function.hpp"
30 #include <set>
31 #include <stack>
32 #include "code_generator.hpp"
33 #include "importer.hpp"
34 #include "options.hpp"
35 #include "shared_object.hpp"
36 #include "timing.hpp"
37 #ifdef CASADI_WITH_THREAD
38 #ifdef CASADI_WITH_THREAD_MINGW
39 #include <mingw.mutex.h>
40 #else // CASADI_WITH_THREAD_MINGW
41 #include <mutex>
42 #endif // CASADI_WITH_THREAD_MINGW
43 #endif //CASADI_WITH_THREAD
44 
45 // This macro is for documentation purposes
46 #define INPUTSCHEME(name)
47 
48 // This macro is for documentation purposes
49 #define OUTPUTSCHEME(name)
50 
52 
53 namespace casadi {
54  template<typename T>
55  std::vector<std::pair<std::string, T>> zip(const std::vector<std::string>& id,
56  const std::vector<T>& mat) {
57  casadi_assert_dev(id.size()==mat.size());
58  std::vector<std::pair<std::string, T>> r(id.size());
59  for (casadi_uint i=0; i<r.size(); ++i) r[i] = std::make_pair(id[i], mat[i]);
60  return r;
61  }
62 
66  struct CASADI_EXPORT ProtoFunctionMemory {
67  // Function specific statistics
68  std::map<std::string, FStats> fstats;
69 
70  // Short-hand for "total" fstats
71  FStats* t_total;
72 
73  // Add a statistic
74  void add_stat(const std::string& s) {
75  bool added = fstats.insert(std::make_pair(s, FStats())).second;
76  casadi_assert(added, "Duplicate stat: '" + s + "'");
77  }
78  };
79 
83  struct CASADI_EXPORT FunctionMemory : public ProtoFunctionMemory {
84  bool stats_available;
85  FunctionMemory() : stats_available(false) {}
86  };
87 
94  class CASADI_EXPORT ProtoFunction : public SharedObjectInternal {
95  public:
99  ProtoFunction(const std::string& name);
100 
104  ~ProtoFunction() override = 0;
105 
111  void construct(const Dict& opts);
112 
114 
117  static const Options options_;
118  virtual const Options& get_options() const { return options_;}
120 
122  virtual Dict generate_options(const std::string& target) const;
123 
127  void print_options(std::ostream &stream) const;
128 
132  void print_option(const std::string &name, std::ostream &stream) const;
133 
137  bool has_option(const std::string &option_name) const;
138 
142  virtual void change_option(const std::string& option_name, const GenericType& option_value);
143 
152  virtual void init(const Dict& opts);
153 
160  virtual void finalize();
161 
163  int checkout() const;
164 
166  void release(int mem) const;
167 
169  void* memory(int ind) const;
170 
172  bool has_memory(int ind) const;
173 
179  virtual void check_mem_count(casadi_int n) const { }
180 
184  virtual void* alloc_mem() const { return new ProtoFunctionMemory(); }
185 
189  virtual int init_mem(void* mem) const;
190 
194  virtual void free_mem(void *mem) const { delete static_cast<ProtoFunctionMemory*>(mem); }
195 
197  virtual Dict get_stats(void* mem) const;
198 
202  void clear_mem();
203 
207  void print(const char* fmt, ...) const;
208 
212  void sprint(char* buf, size_t buf_sz, const char* fmt, ...) const;
213 
217  void format_time(char* buffer, double time) const;
218 
222  void print_time(const std::map<std::string, FStats>& fstats) const;
223 
227  void serialize(SerializingStream &s) const;
228 
232  virtual void serialize_body(SerializingStream &s) const;
236  virtual void serialize_type(SerializingStream &s) const {}
237 
241  virtual std::string serialize_base_function() const {
242  return class_name();
243  }
244 
246  std::string name_;
247 
249  bool verbose_;
250 
251  // Print timing statistics
252  bool print_time_;
253 
254  // Print timing statistics
255  bool record_time_;
256 
258  bool regularity_check_;
259 
261  bool error_on_fail_;
262 
263  protected:
267  explicit ProtoFunction(DeserializingStream& s);
268 
269 #ifdef CASADI_WITH_THREAD
271  mutable std::mutex mtx_;
272 #endif // CASADI_WITH_THREAD
273 
274  private:
276  mutable std::vector<void*> mem_;
277 
279  mutable std::stack<int> unused_;
280  };
281 
288  class CASADI_EXPORT FunctionInternal : public ProtoFunction {
289  friend class Function;
290  public:
294  FunctionInternal(const std::string& name);
295 
299  ~FunctionInternal() override = 0;
300 
304  virtual std::string getAdaptorSolverName() const { return ""; }
305 
307 
310  static const Options options_;
311  const Options& get_options() const override { return options_;}
313 
315  Dict generate_options(const std::string& target) const override;
316 
320  void change_option(const std::string& option_name, const GenericType& option_value) override;
321 
325  void reset_dump_count();
326 
330  void init(const Dict& opts) override;
331 
335  void finalize() override;
336 
340  void* alloc_mem() const override { return new FunctionMemory(); }
341 
345  void free_mem(void *mem) const override { delete static_cast<FunctionMemory*>(mem); }
346 
348  Dict get_stats(void* mem) const override;
349 
353  Function self() const { return shared_from_this<Function>();}
354 
355  // Factory
356  virtual Function factory(const std::string& name,
357  const std::vector<std::string>& s_in,
358  const std::vector<std::string>& s_out,
359  const Function::AuxOut& aux,
360  const Dict& opts) const;
361 
362  // Get list of dependency functions
363  virtual std::vector<std::string> get_function() const;
364 
365  // Get a dependency function
366  virtual const Function& get_function(const std::string &name) const;
367 
368  // Check if a particular dependency exists
369  virtual bool has_function(const std::string& fname) const {return false;}
370 
371  // Add embedded function to map, helper function
372  void add_embedded(std::map<FunctionInternal*, Function>& all_fun,
373  const Function& dep, casadi_int max_depth) const;
374 
375  // Get all embedded functions, recursively
376  virtual void find(std::map<FunctionInternal*, Function>& all_fun, casadi_int max_depth) const {}
377 
386  virtual std::vector<bool> which_depends(const std::string& s_in,
387  const std::vector<std::string>& s_out,
388  casadi_int order, bool tr=false) const;
389 
391 
394  virtual bool has_spfwd() const { return false;}
395  virtual bool has_sprev() const { return false;}
397 
399 
402  int eval_gen(const double** arg, double** res, casadi_int* iw, double* w, void* mem,
403  bool always_inline, bool never_inline) const;
404  virtual int eval(const double** arg, double** res, casadi_int* iw, double* w, void* mem) const;
406 
410  virtual int eval_sx(const SXElem** arg, SXElem** res,
411  casadi_int* iw, SXElem* w, void* mem, bool always_inline, bool never_inline) const;
412 
416  virtual void eval_mx(const MXVector& arg, MXVector& res,
417  bool always_inline, bool never_inline) const;
418 
420 
423  virtual std::vector<DM> eval_dm(const std::vector<DM>& arg) const;
424  virtual bool has_eval_dm() const { return false;}
426 
428 
431  int eval_gen(const SXElem** arg, SXElem** res, casadi_int* iw, SXElem* w, void* mem,
432  bool always_inline, bool never_inline) const {
433  return eval_sx(arg, res, iw, w, mem, always_inline, never_inline);
434  }
435  int eval_gen(const bvec_t** arg, bvec_t** res, casadi_int* iw, bvec_t* w, void* mem,
436  bool always_inline, bool never_inline) const {
437  return sp_forward(arg, res, iw, w, mem);
438  }
440 
442 
445  void call_gen(const MXVector& arg, MXVector& res, casadi_int npar,
446  bool always_inline, bool never_inline) const;
447 
448  template<typename D>
449  void call_gen(const std::vector<Matrix<D> >& arg, std::vector<Matrix<D> >& res,
450  casadi_int npar, bool always_inline, bool never_inline) const;
452 
456  template<typename M>
457  void call(const std::vector<M>& arg, std::vector<M>& res,
458  bool always_inline, bool never_inline) const;
459 
461 
466  static bool check_mat(const Sparsity& arg, const Sparsity& inp, casadi_int& npar);
468 
470 
478  template<typename M>
479  void check_arg(const std::vector<M>& arg, casadi_int& npar) const;
481 
483 
491  template<typename M>
492  void check_res(const std::vector<M>& res, casadi_int& npar) const;
494 
503  template<typename M> bool
504  matching_arg(const std::vector<M>& arg, casadi_int& npar) const;
505 
514  template<typename M> bool
515  matching_res(const std::vector<M>& arg, casadi_int& npar) const;
516 
520  template<typename M> std::vector<M>
521  replace_arg(const std::vector<M>& arg, casadi_int npar) const;
522 
526  template<typename M> std::vector<M>
527  project_arg(const std::vector<M>& arg, casadi_int npar) const;
528 
532  template<typename M> std::vector<M>
533  project_res(const std::vector<M>& arg, casadi_int npar) const;
534 
538  template<typename M> std::vector<M>
539  replace_res(const std::vector<M>& res, casadi_int npar) const;
540 
544  template<typename M> std::vector<std::vector<M>>
545  replace_fseed(const std::vector<std::vector<M>>& fseed, casadi_int npar) const;
546 
550  template<typename M> std::vector<std::vector<M>>
551  replace_aseed(const std::vector<std::vector<M>>& aseed, casadi_int npar) const;
552 
557  template<typename M>
558  std::map<std::string, M> convert_arg(const std::vector<M>& arg) const;
559  template<typename M>
560  std::vector<M> convert_arg(const std::map<std::string, M>& arg) const;
561  template<typename M>
562  std::map<std::string, M> convert_res(const std::vector<M>& res) const;
563  template<typename M>
564  std::vector<M> convert_res(const std::map<std::string, M>& res) const;
566 
571  std::vector<double> nz_in(const std::vector<DM>& arg) const;
572  std::vector<double> nz_out(const std::vector<DM>& res) const;
573  std::vector<DM> nz_in(const std::vector<double>& arg) const;
574  std::vector<DM> nz_out(const std::vector<double>& res) const;
576 
578 
581  virtual void call_forward(const std::vector<MX>& arg, const std::vector<MX>& res,
582  const std::vector<std::vector<MX> >& fseed,
583  std::vector<std::vector<MX> >& fsens,
584  bool always_inline, bool never_inline) const;
585  virtual void call_forward(const std::vector<SX>& arg, const std::vector<SX>& res,
586  const std::vector<std::vector<SX> >& fseed,
587  std::vector<std::vector<SX> >& fsens,
588  bool always_inline, bool never_inline) const;
590 
592 
595  virtual void call_reverse(const std::vector<MX>& arg, const std::vector<MX>& res,
596  const std::vector<std::vector<MX> >& aseed,
597  std::vector<std::vector<MX> >& asens,
598  bool always_inline, bool never_inline) const;
599  virtual void call_reverse(const std::vector<SX>& arg, const std::vector<SX>& res,
600  const std::vector<std::vector<SX> >& aseed,
601  std::vector<std::vector<SX> >& asens,
602  bool always_inline, bool never_inline) const;
604 
608  std::vector<MX> mapsum_mx(const std::vector<MX > &arg, const std::string& parallelization);
609 
613  virtual bool uses_output() const {return false;}
614 
616 
619  Function jacobian() const;
620  virtual bool has_jacobian() const { return false;}
621  virtual Function get_jacobian(const std::string& name,
622  const std::vector<std::string>& inames,
623  const std::vector<std::string>& onames,
624  const Dict& opts) const;
626 
628 
632  Sparsity& jac_sparsity(casadi_int oind, casadi_int iind, bool compact, bool symmetric) const;
633  virtual bool has_jac_sparsity(casadi_int oind, casadi_int iind) const { return false;}
634  virtual Sparsity get_jac_sparsity(casadi_int oind, casadi_int iind, bool symmetric) const;
636 
638  static std::string forward_name(const std::string& fcn, casadi_int nfwd) {
639  return "fwd" + str(nfwd) + "_" + fcn;
640  }
641 
643  std::string diff_prefix(const std::string& prefix) const;
644 
646 
653  Function forward(casadi_int nfwd) const;
654  virtual bool has_forward(casadi_int nfwd) const { return false;}
655  virtual Function get_forward(casadi_int nfwd, const std::string& name,
656  const std::vector<std::string>& inames,
657  const std::vector<std::string>& onames,
658  const Dict& opts) const;
660 
662  static std::string reverse_name(const std::string& fcn, casadi_int nadj) {
663  return "adj" + str(nadj) + "_" + fcn;
664  }
665 
667 
674  Function reverse(casadi_int nadj) const;
675  virtual bool has_reverse(casadi_int nadj) const { return false;}
676  virtual Function get_reverse(casadi_int nadj, const std::string& name,
677  const std::vector<std::string>& inames,
678  const std::vector<std::string>& onames,
679  const Dict& opts) const;
681 
685  template<typename MatType>
686  static MatType ensure_stacked(const MatType& v, const Sparsity& sp, casadi_int n);
687 
691  virtual Function slice(const std::string& name, const std::vector<casadi_int>& order_in,
692  const std::vector<casadi_int>& order_out, const Dict& opts) const;
693 
697  virtual const Function& oracle() const;
698 
702  bool has_derivative() const;
703 
707  virtual double ad_weight() const;
708 
714  virtual double sp_weight() const;
715 
717 
720  virtual const SX sx_in(casadi_int ind) const;
721  virtual const SX sx_out(casadi_int ind) const;
722  virtual const std::vector<SX> sx_in() const;
723  virtual const std::vector<SX> sx_out() const;
724  virtual const MX mx_in(casadi_int ind) const;
725  virtual const MX mx_out(casadi_int ind) const;
726  virtual const std::vector<MX> mx_in() const;
727  virtual const std::vector<MX> mx_out() const;
728  const DM dm_in(casadi_int ind) const;
729  const DM dm_out(casadi_int ind) const;
730  const std::vector<DM> dm_in() const;
731  const std::vector<DM> dm_out() const;
733 
735  virtual std::vector<MX> free_mx() const;
736 
738  virtual std::vector<SX> free_sx() const;
739 
743  virtual bool has_free() const { return false;}
744 
748  virtual void generate_lifted(Function& vdef_fcn, Function& vinit_fcn) const;
749 
753  virtual casadi_int n_instructions() const;
754 
758  virtual casadi_int instruction_id(casadi_int k) const;
759 
763  virtual std::vector<casadi_int> instruction_input(casadi_int k) const;
764 
768  virtual double instruction_constant(casadi_int k) const;
769 
773  virtual std::vector<casadi_int> instruction_output(casadi_int k) const;
774 
778  virtual casadi_int n_nodes() const;
779 
783  virtual MX instruction_MX(casadi_int k) const;
784 
788  virtual SX instructions_sx() const;
789 
793  Function wrap() const;
794 
798  Function wrap_as_needed(const Dict& opts) const;
799 
803  Dict cache() const;
804 
808  bool incache(const std::string& fname, Function& f, const std::string& suffix="") const;
809 
813  void tocache(const Function& f, const std::string& suffix="") const;
814 
815 
819  void tocache_if_missing(Function& f, const std::string& suffix="") const;
820 
824  void codegen(CodeGenerator& g, const std::string& fname) const;
825 
829  void codegen_meta(CodeGenerator& g) const;
830 
834  void codegen_sparsities(CodeGenerator& g) const;
835 
839  virtual std::string codegen_name(const CodeGenerator& g, bool ns=true) const;
840 
844  std::string codegen_mem(CodeGenerator& g, const std::string& index="mem") const;
845 
849  virtual void codegen_incref(CodeGenerator& g) const {}
850 
854  virtual void codegen_decref(CodeGenerator& g) const {}
855 
859  virtual void codegen_alloc_mem(CodeGenerator& g) const;
860 
864  virtual void codegen_init_mem(CodeGenerator& g) const;
865 
869  virtual void codegen_free_mem(CodeGenerator& g) const {}
870 
874  virtual void codegen_checkout(CodeGenerator& g) const;
875 
879  virtual void codegen_release(CodeGenerator& g) const;
880 
884  std::string signature(const std::string& fname) const;
885 
889  std::string signature_unrolled(const std::string& fname) const;
890 
894  virtual void codegen_declarations(CodeGenerator& g) const;
895 
899  virtual void codegen_body(CodeGenerator& g) const;
900 
904  virtual std::string codegen_mem_type() const { return ""; }
905 
909  virtual std::string generate_dependencies(const std::string& fname, const Dict& opts) const;
910 
914  virtual bool has_codegen() const { return false;}
915 
919  virtual void jit_dependencies(const std::string& fname) {}
920 
924  virtual void export_code(const std::string& lang,
925  std::ostream &stream, const Dict& options) const;
926 
930  void serialize_type(SerializingStream &s) const override;
931 
935  void serialize_body(SerializingStream &s) const override;
936 
940  void disp(std::ostream& stream, bool more) const override;
941 
945  virtual void disp_more(std::ostream& stream) const {}
946 
950  std::string definition() const;
951 
955  void print_dimensions(std::ostream &stream) const;
956 
960  virtual std::vector<std::string> get_free() const;
961 
965  void get_partition(casadi_int iind, casadi_int oind, Sparsity& D1, Sparsity& D2,
966  bool compact, bool symmetric,
967  bool allow_forward, bool allow_reverse) const;
968 
970 
973  casadi_int nnz_in() const;
974  casadi_int nnz_in(casadi_int ind) const { return sparsity_in(ind).nnz(); }
975  casadi_int nnz_out() const;
976  casadi_int nnz_out(casadi_int ind) const { return sparsity_out(ind).nnz(); }
978 
980 
983  casadi_int numel_in() const;
984  casadi_int numel_in(casadi_int ind) const { return sparsity_in(ind).numel(); }
985  casadi_int numel_out(casadi_int ind) const { return sparsity_out(ind).numel(); }
986  casadi_int numel_out() const;
988 
990 
993  casadi_int size1_in(casadi_int ind) const { return sparsity_in(ind).size1(); }
994  casadi_int size2_in(casadi_int ind) const { return sparsity_in(ind).size2(); }
995  casadi_int size1_out(casadi_int ind) const { return sparsity_out(ind).size1(); }
996  casadi_int size2_out(casadi_int ind) const { return sparsity_out(ind).size2(); }
997  std::pair<casadi_int, casadi_int> size_in(casadi_int ind) const {
998  return sparsity_in(ind).size();
999  }
1000  std::pair<casadi_int, casadi_int> size_out(casadi_int ind) const {
1001  return sparsity_out(ind).size();
1002  }
1004 
1006 
1009  const Sparsity& sparsity_in(casadi_int ind) const { return sparsity_in_.at(ind); }
1010  const Sparsity& sparsity_out(casadi_int ind) const { return sparsity_out_.at(ind); }
1012 
1014 
1017  bool all_scalar() const;
1018 
1020  virtual bool jac_is_symm(casadi_int oind, casadi_int iind) const;
1021 
1023  Sparsity to_compact(casadi_int oind, casadi_int iind, const Sparsity& sp) const;
1024 
1026  Sparsity from_compact(casadi_int oind, casadi_int iind, const Sparsity& sp) const;
1027 
1029  template<bool fwd>
1030  Sparsity get_jac_sparsity_gen(casadi_int oind, casadi_int iind) const;
1031 
1033  Sparsity get_jac_sparsity_hierarchical(casadi_int oind, casadi_int iind) const;
1034 
1038  Sparsity get_jac_sparsity_hierarchical_symm(casadi_int oind, casadi_int iind) const;
1039 
1041  virtual std::vector<MX> symbolic_output(const std::vector<MX>& arg) const;
1042 
1044 
1047  virtual size_t get_n_in();
1048  virtual size_t get_n_out();
1050 
1051 
1053 
1056  virtual std::string get_name_in(casadi_int i);
1057  virtual std::string get_name_out(casadi_int i);
1059 
1063  virtual double get_default_in(casadi_int ind) const {
1064  return 0;
1065  }
1066 
1070  virtual double get_max_in(casadi_int ind) const {
1071  return inf;
1072  }
1073 
1077  virtual double get_min_in(casadi_int ind) const {
1078  return -inf;
1079  }
1080 
1081  virtual std::vector<double> get_nominal_in(casadi_int ind) const {
1082  return std::vector<double>(nnz_in(ind), 1.);
1083  }
1084 
1085  virtual std::vector<double> get_nominal_out(casadi_int ind) const {
1086  return std::vector<double>(nnz_out(ind), 1.);
1087  }
1088 
1092  virtual double get_reltol() const {
1093  return eps;
1094  }
1095 
1099  virtual double get_abstol() const {
1100  return eps;
1101  }
1102 
1106  virtual Sparsity get_sparsity_in(casadi_int i);
1107 
1111  virtual Sparsity get_sparsity_out(casadi_int i);
1112 
1116  virtual bool get_diff_in(casadi_int i) { return true; }
1117 
1121  virtual bool get_diff_out(casadi_int i) { return true; }
1122 
1126  casadi_int index_in(const std::string &name) const {
1127  for (casadi_int i=0; i<name_in_.size(); ++i) {
1128  if (name_in_[i]==name) return i;
1129  }
1130  casadi_error("FunctionInternal::index_in: could not find entry \""
1131  + name + "\". Available names are: " + str(name_in_) + ".");
1132  return -1;
1133  }
1134 
1138  casadi_int index_out(const std::string &name) const {
1139  for (casadi_int i=0; i<name_out_.size(); ++i) {
1140  if (name_out_[i]==name) return i;
1141  }
1142  casadi_error("FunctionInternal::index_out: could not find entry \""
1143  + name + "\". Available names are: " + str(name_out_) + ".");
1144  return -1;
1145  }
1146 
1150  virtual int sp_forward(const bvec_t** arg, bvec_t** res,
1151  casadi_int* iw, bvec_t* w, void* mem) const;
1152 
1156  virtual int sp_forward_block(const bvec_t** arg, bvec_t** res,
1157  casadi_int* iw, bvec_t* w, void* mem, casadi_int oind, casadi_int iind) const;
1158 
1162  virtual int sp_reverse(bvec_t** arg, bvec_t** res, casadi_int* iw, bvec_t* w, void* mem) const;
1163 
1167  void sz_work(size_t& sz_arg, size_t& sz_res, size_t& sz_iw, size_t& sz_w) const;
1168 
1172  size_t sz_arg() const { return sz_arg_per_ + sz_arg_tmp_;}
1173 
1177  size_t sz_res() const { return sz_res_per_ + sz_res_tmp_;}
1178 
1182  size_t sz_iw() const { return sz_iw_per_ + sz_iw_tmp_;}
1183 
1187  size_t sz_w() const { return sz_w_per_ + sz_w_tmp_;}
1188 
1193  virtual size_t codegen_sz_arg(const CodeGenerator& g) const;
1194  virtual size_t codegen_sz_res(const CodeGenerator& g) const;
1195  virtual size_t codegen_sz_iw(const CodeGenerator& g) const;
1196  virtual size_t codegen_sz_w(const CodeGenerator& g) const;
1198 
1202  void alloc_arg(size_t sz_arg, bool persistent=false);
1203 
1207  void alloc_res(size_t sz_res, bool persistent=false);
1208 
1212  void alloc_iw(size_t sz_iw, bool persistent=false);
1213 
1217  void alloc_w(size_t sz_w, bool persistent=false);
1218 
1222  void alloc(const Function& f, bool persistent=false, int num_threads=1);
1223 
1227  virtual void set_work(void* mem, const double**& arg, double**& res,
1228  casadi_int*& iw, double*& w) const {}
1229 
1233  virtual void set_temp(void* mem, const double** arg, double** res,
1234  casadi_int* iw, double* w) const {}
1235 
1239  void setup(void* mem, const double** arg, double** res, casadi_int* iw, double* w) const;
1240 
1242 
1245  virtual bool fwdViaJac(casadi_int nfwd) const;
1246  virtual bool adjViaJac(casadi_int nadj) const;
1248 
1250  virtual Dict info() const;
1251 
1255  Function map(casadi_int n, const std::string& parallelization) const;
1256 
1260  void generate_in(const std::string& fname, const double** arg) const;
1261  void generate_out(const std::string& fname, double** res) const;
1262 
1263  bool always_inline_, never_inline_;
1264 
1266  size_t n_in_, n_out_;
1267 
1269  std::vector<bool> is_diff_in_, is_diff_out_;
1270 
1272  std::vector<Sparsity> sparsity_in_, sparsity_out_;
1273 
1275  std::vector<std::string> name_in_, name_out_;
1276 
1280  bool jit_;
1281 
1285  bool jit_cleanup_;
1286 
1290  std::string jit_serialize_;
1291 
1295  std::string jit_name_;
1296 
1297  std::string jit_base_name_;
1298 
1302  bool jit_temp_suffix_;
1303 
1307  eval_t eval_;
1308 
1312  casadi_checkout_t checkout_;
1313 
1317  casadi_release_t release_;
1318 
1322  Dict stats_;
1323 
1327  bool has_refcount_;
1328 
1332  Dict cache_init_;
1333 
1335  mutable WeakCache<std::string, Function> cache_;
1336 
1338  mutable std::vector<Sparsity> jac_sparsity_[2];
1339 
1340 #ifdef CASADI_WITH_THREADSAFE_SYMBOLICS
1342  mutable std::mutex jac_sparsity_mtx_;
1343 #endif // CASADI_WITH_THREADSAFE_SYMBOLICS
1344 
1346  Function derivative_of_;
1347 
1349  void* user_data_;
1350 
1352  std::string compiler_plugin_;
1353  Importer compiler_;
1354  Dict jit_options_;
1355 
1357  double jac_penalty_;
1358 
1359  // Types of derivative calculation permitted
1360  bool enable_forward_, enable_reverse_, enable_jacobian_, enable_fd_;
1361  bool enable_forward_op_, enable_reverse_op_, enable_jacobian_op_, enable_fd_op_;
1362 
1364  double ad_weight_, ad_weight_sp_;
1365 
1367  casadi_int max_num_dir_;
1368 
1370  bool inputs_check_;
1371 
1372  // Finite difference step
1373  Dict fd_options_;
1374 
1375  // Finite difference step size
1376  double fd_step_;
1377 
1378  // Finite difference method
1379  std::string fd_method_;
1380 
1381  // Print input/output
1382  bool print_in_;
1383  bool print_out_;
1384 
1385  // Print canonical form
1386  bool print_canonical_;
1387 
1388  // Warn when number of inputs or outputs exceed this value
1389  casadi_int max_io_;
1390 
1391  // Dump input/output
1392  bool dump_in_, dump_out_, dump_;
1393 
1394  // Directory to dump to
1395  std::string dump_dir_;
1396 
1397  // Format to dump with
1398  std::string dump_format_;
1399 
1400  // Forward/reverse/Jacobian options
1401  Dict forward_options_, reverse_options_, jacobian_options_, der_options_;
1402 
1403  // Store a reference to a custom Jacobian
1404  Function custom_jacobian_;
1405 
1406  // Counter for unique names for dumping inputs and output
1407 #ifdef CASADI_WITH_THREAD
1408  mutable std::atomic<casadi_int> dump_count_;
1409 #else
1410  mutable casadi_int dump_count_;
1411 #endif // CASADI_WITH_THREAD
1412 
1416  virtual bool is_a(const std::string& type, bool recursive) const;
1417 
1421  virtual void merge(const std::vector<MX>& arg,
1422  std::vector<MX>& subs_from, std::vector<MX>& subs_to) const;
1423 
1427  template<typename MatType>
1428  static bool purgable(const std::vector<MatType>& seed);
1429 
1433  template<typename MatType>
1434  std::vector<std::vector<MatType> >
1435  fwd_seed(casadi_int nfwd) const;
1436 
1440  template<typename MatType>
1441  std::vector<std::vector<MatType> >
1442  symbolicAdjSeed(casadi_int nadj, const std::vector<MatType>& v) const;
1443 
1444  static std::string string_from_UnifiedReturnStatus(UnifiedReturnStatus status);
1445 
1449  explicit FunctionInternal(DeserializingStream& e);
1450 
1454  static Function deserialize(DeserializingStream& s);
1455  static std::map<std::string, ProtoFunction* (*)(DeserializingStream&)> deserialize_map;
1456 
1460  void print_in(std::ostream &stream, const double** arg, bool truncate) const;
1461 
1465  void print_out(std::ostream &stream, double** res, bool truncate) const;
1466 
1470  static void print_canonical(std::ostream &stream, const Sparsity& sp, const double* nz);
1471 
1475  static void print_canonical(std::ostream &stream, casadi_int sz, const double* nz);
1476 
1480  static void print_canonical(std::ostream &stream, double a);
1481 
1482  protected:
1486  void set_jac_sparsity(casadi_int oind, casadi_int iind, const Sparsity& sp);
1487 
1488  private:
1489  // @{
1491  casadi_int get_dump_id() const;
1492  void dump_in(casadi_int id, const double** arg) const;
1493  void dump_out(casadi_int id, double** res) const;
1494  void dump() const;
1495  // @}
1496 
1500  size_t sz_arg_per_, sz_res_per_, sz_iw_per_, sz_w_per_;
1501 
1505  size_t sz_arg_tmp_, sz_res_tmp_, sz_iw_tmp_, sz_w_tmp_;
1506  };
1507 
1508  // Template implementations
1509  template<typename MatType>
1510  bool FunctionInternal::purgable(const std::vector<MatType>& v) {
1511  for (auto i=v.begin(); i!=v.end(); ++i) {
1512  if (!i->is_zero()) return false;
1513  }
1514  return true;
1515  }
1516 
1517  template<typename MatType>
1518  std::vector<std::vector<MatType> >
1519  FunctionInternal::
1520  fwd_seed(casadi_int nfwd) const {
1521  std::vector<std::vector<MatType>> fseed(nfwd);
1522  for (casadi_int dir=0; dir<nfwd; ++dir) {
1523  fseed[dir].resize(n_in_);
1524  for (casadi_int iind=0; iind<n_in_; ++iind) {
1525  std::string n = "f" + str(dir) + "_" + name_in_[iind];
1526  Sparsity sp = is_diff_in_[iind] ? sparsity_in(iind) : Sparsity(size_in(iind));
1527  fseed[dir][iind] = MatType::sym(n, sp);
1528  }
1529  }
1530  return fseed;
1531  }
1532 
1533  template<typename MatType>
1534  std::vector<std::vector<MatType> >
1535  FunctionInternal::
1536  symbolicAdjSeed(casadi_int nadj, const std::vector<MatType>& v) const {
1537  std::vector<std::vector<MatType> > aseed(nadj, v);
1538  for (casadi_int dir=0; dir<nadj; ++dir) {
1539  // Replace symbolic inputs
1540  casadi_int oind=0;
1541  for (typename std::vector<MatType>::iterator i=aseed[dir].begin();
1542  i!=aseed[dir].end();
1543  ++i, ++oind) {
1544  // Name of the adjoint seed
1545  std::stringstream ss;
1546  ss << "a";
1547  if (nadj>1) ss << dir << "_";
1548  ss << oind;
1549 
1550  // Save to matrix
1551  *i = MatType::sym(ss.str(), is_diff_out_[oind] ? i->sparsity() : Sparsity(i->size()));
1552 
1553  }
1554  }
1555  return aseed;
1556  }
1557 
1558  template<typename M>
1559  void FunctionInternal::call(const std::vector<M>& arg, std::vector<M>& res,
1560  bool always_inline, bool never_inline) const {
1561  // If all inputs are scalar ...
1562  if (all_scalar()) {
1563  // ... and some arguments are matrix-valued with matching dimensions ...
1564  bool matrix_call = false;
1565  std::pair<casadi_int, casadi_int> sz;
1566  for (auto&& a : arg) {
1567  if (!a.is_scalar() && !a.is_empty()) {
1568  if (!matrix_call) {
1569  // Matrix call
1570  matrix_call = true;
1571  sz = a.size();
1572  } else if (a.size()!=sz) {
1573  // Not same dimensions
1574  matrix_call = false;
1575  break;
1576  }
1577  }
1578  }
1579 
1580  // ... then, call multiple times
1581  if (matrix_call) {
1582  // Start with zeros
1583  res.resize(n_out_);
1584  M z = M::zeros(sz);
1585  for (auto&& a : res) a = z;
1586  // Call multiple times
1587  std::vector<M> arg1 = arg, res1;
1588  for (casadi_int c=0; c<sz.second; ++c) {
1589  for (casadi_int r=0; r<sz.first; ++r) {
1590  // Get scalar arguments
1591  for (casadi_int i=0; i<arg.size(); ++i) {
1592  if (arg[i].size()==sz) arg1[i] = arg[i](r, c);
1593  }
1594  // Call recursively with scalar arguments
1595  call(arg1, res1, always_inline, never_inline);
1596  // Get results
1597  casadi_assert_dev(res.size() == res1.size());
1598  for (casadi_int i=0; i<res.size(); ++i) res[i](r, c) = res1[i];
1599  }
1600  }
1601  // All elements assigned
1602  return;
1603  }
1604  }
1605 
1606  // Check if inputs need to be replaced
1607  casadi_int npar = 1;
1608  if (!matching_arg(arg, npar)) {
1609  return call(replace_arg(arg, npar), res, always_inline, never_inline);
1610  }
1611 
1612  // Call the type-specific method
1613  call_gen(arg, res, npar, always_inline, never_inline);
1614  }
1615 
1616  template<typename M>
1617  std::vector<M> FunctionInternal::
1618  project_arg(const std::vector<M>& arg, casadi_int npar) const {
1619  casadi_assert_dev(arg.size()==n_in_);
1620 
1621  // Which arguments require mapped evaluation
1622  std::vector<bool> mapped(n_in_);
1623  for (casadi_int i=0; i<n_in_; ++i) {
1624  mapped[i] = arg[i].size2()!=size2_in(i);
1625  }
1626 
1627  // Check if matching input sparsity
1628  std::vector<bool> matching(n_in_);
1629  bool any_mismatch = false;
1630  for (casadi_int i=0; i<n_in_; ++i) {
1631  if (mapped[i]) {
1632  matching[i] = arg[i].sparsity().is_stacked(sparsity_in(i), npar);
1633  } else {
1634  matching[i] = arg[i].sparsity()==sparsity_in(i);
1635  }
1636  any_mismatch = any_mismatch || !matching[i];
1637  }
1638 
1639  // Correct input sparsity
1640  if (any_mismatch) {
1641  std::vector<M> arg2(arg);
1642  for (casadi_int i=0; i<n_in_; ++i) {
1643  if (!matching[i]) {
1644  if (mapped[i]) {
1645  arg2[i] = project(arg2[i], repmat(sparsity_in(i), 1, npar));
1646  } else {
1647  arg2[i] = project(arg2[i], sparsity_in(i));
1648  }
1649  }
1650  }
1651  return arg2;
1652  }
1653  return arg;
1654  }
1655 
1656  template<typename M>
1657  std::vector<M> FunctionInternal::
1658  project_res(const std::vector<M>& arg, casadi_int npar) const {
1659  return arg;
1660  }
1661 
1662  template<typename D>
1663  void FunctionInternal::
1664  call_gen(const std::vector<Matrix<D> >& arg, std::vector<Matrix<D> >& res,
1665  casadi_int npar, bool always_inline, bool never_inline) const {
1666  std::vector< Matrix<D> > arg2 = project_arg(arg, npar);
1667 
1668  // Which arguments require mapped evaluation
1669  std::vector<bool> mapped(n_in_);
1670  for (casadi_int i=0; i<n_in_; ++i) {
1671  mapped[i] = arg[i].size2()!=size2_in(i);
1672  }
1673 
1674  // Allocate results
1675  res.resize(n_out_);
1676  for (casadi_int i=0; i<n_out_; ++i) {
1677  if (!res[i].sparsity().is_stacked(sparsity_out(i), npar)) {
1678  res[i] = Matrix<D>::zeros(repmat(sparsity_out(i), 1, npar));
1679  }
1680  }
1681 
1682  // Allocate temporary memory if needed
1683  std::vector<casadi_int> iw_tmp(sz_iw());
1684  std::vector<D> w_tmp(sz_w());
1685 
1686  // Get pointers to input arguments
1687  std::vector<const D*> argp(sz_arg());
1688  for (casadi_int i=0; i<n_in_; ++i) argp[i]=get_ptr(arg2[i]);
1689 
1690  // Get pointers to output arguments
1691  std::vector<D*> resp(sz_res());
1692  for (casadi_int i=0; i<n_out_; ++i) resp[i]=get_ptr(res[i]);
1693 
1694  // For all parallel calls
1695  for (casadi_int p=0; p<npar; ++p) {
1696  // Call memory-less
1697  if (eval_gen(get_ptr(argp), get_ptr(resp),
1698  get_ptr(iw_tmp), get_ptr(w_tmp), memory(0),
1699  always_inline, never_inline)) {
1700  if (error_on_fail_) casadi_error("Evaluation failed");
1701  }
1702  // Update offsets
1703  if (p==npar-1) break;
1704  for (casadi_int i=0; i<n_in_; ++i) if (mapped[i]) argp[i] += nnz_in(i);
1705  for (casadi_int i=0; i<n_out_; ++i) resp[i] += nnz_out(i);
1706  }
1707  }
1708 
1709  template<typename M>
1710  void FunctionInternal::check_arg(const std::vector<M>& arg, casadi_int& npar) const {
1711  casadi_assert(arg.size()==n_in_, "Incorrect number of inputs: Expected "
1712  + str(n_in_) + ", got " + str(arg.size()));
1713  for (casadi_int i=0; i<n_in_; ++i) {
1714  if (!check_mat(arg[i].sparsity(), sparsity_in(i), npar)) {
1715  // Dimensions
1716  std::string d_arg = str(arg[i].size1()) + "-by-" + str(arg[i].size2());
1717  std::string d_in = str(size1_in(i)) + "-by-" + str(size2_in(i));
1718  std::string e = "Input " + str(i) + " (" + name_in_[i] + ") has mismatching shape. "
1719  "Got " + d_arg + ". Allowed dimensions, in general, are:\n"
1720  " - The input dimension N-by-M (here " + d_in + ")\n"
1721  " - A scalar, i.e. 1-by-1\n"
1722  " - M-by-N if N=1 or M=1 (i.e. a transposed vector)\n"
1723  " - N-by-M1 if K*M1=M for some K (argument repeated horizontally)\n";
1724  if (npar!=-1) {
1725  e += " - N-by-P*M, indicating evaluation with multiple arguments (P must be a "
1726  "multiple of " + str(npar) + " for consistency with previous inputs)";
1727  }
1728  casadi_error(e);
1729  }
1730  }
1731  }
1732 
1733  template<typename M>
1734  void FunctionInternal::check_res(const std::vector<M>& res, casadi_int& npar) const {
1735  casadi_assert(res.size()==n_out_, "Incorrect number of outputs: Expected "
1736  + str(n_out_) + ", got " + str(res.size()));
1737  for (casadi_int i=0; i<n_out_; ++i) {
1738  casadi_assert(check_mat(res[i].sparsity(), sparsity_out(i), npar),
1739  "Output " + str(i) + " (" + name_out_[i] + ") has mismatching shape. "
1740  "Expected " + str(size_out(i)) + ", got " + str(res[i].size()));
1741  }
1742  }
1743 
1744  template<typename M>
1745  bool FunctionInternal::matching_arg(const std::vector<M>& arg, casadi_int& npar) const {
1746  check_arg(arg, npar);
1747  for (casadi_int i=0; i<n_in_; ++i) {
1748  if (arg.at(i).size1()!=size1_in(i)) return false;
1749  if (arg.at(i).size2()!=size2_in(i) && arg.at(i).size2()!=npar*size2_in(i)) return false;
1750  }
1751  return true;
1752  }
1753 
1754  template<typename M>
1755  bool FunctionInternal::matching_res(const std::vector<M>& res, casadi_int& npar) const {
1756  check_res(res, npar);
1757  for (casadi_int i=0; i<n_out_; ++i) {
1758  if (res.at(i).size1()!=size1_out(i)) return false;
1759  if (res.at(i).size2()!=size2_out(i) && res.at(i).size2()!=npar*size2_out(i)) return false;
1760  }
1761  return true;
1762  }
1763 
1764  template<typename M>
1765  M replace_mat(const M& arg, const Sparsity& inp, casadi_int npar) {
1766  if (arg.size()==inp.size()) {
1767  // Matching dimensions already
1768  return arg;
1769  } else if (arg.is_empty()) {
1770  // Empty matrix means set zero
1771  return M(inp.size());
1772  } else if (arg.is_scalar()) {
1773  // Scalar assign means set all
1774  return M(inp, arg);
1775  } else if (arg.is_vector() && inp.size()==std::make_pair(arg.size2(), arg.size1())) {
1776  // Transpose vector
1777  return arg.T();
1778  } else if (arg.size1()==inp.size1() && arg.size2()>0 && inp.size2()>0
1779  && inp.size2()%arg.size2()==0) {
1780  // Horizontal repmat
1781  return repmat(arg, 1, inp.size2()/arg.size2());
1782  } else {
1783  casadi_assert_dev(npar!=-1);
1784  // Multiple evaluation
1785  return repmat(arg, 1, (npar*inp.size2())/arg.size2());
1786  }
1787  }
1788 
1789  template<typename M>
1790  std::vector<M> FunctionInternal::
1791  replace_arg(const std::vector<M>& arg, casadi_int npar) const {
1792  std::vector<M> r(arg.size());
1793  for (casadi_int i=0; i<r.size(); ++i) r[i] = replace_mat(arg[i], sparsity_in(i), npar);
1794  return r;
1795  }
1796 
1797  template<typename M>
1798  std::vector<M> FunctionInternal::
1799  replace_res(const std::vector<M>& res, casadi_int npar) const {
1800  std::vector<M> r(res.size());
1801  for (casadi_int i=0; i<r.size(); ++i) r[i] = replace_mat(res[i], sparsity_out(i), npar);
1802  return r;
1803  }
1804 
1805  template<typename M>
1806  std::vector<std::vector<M> > FunctionInternal::
1807  replace_fseed(const std::vector<std::vector<M> >& fseed, casadi_int npar) const {
1808  std::vector<std::vector<M> > r(fseed.size());
1809  for (casadi_int d=0; d<r.size(); ++d) r[d] = replace_arg(fseed[d], npar);
1810  return r;
1811  }
1812 
1813  template<typename M>
1814  std::vector<std::vector<M> > FunctionInternal::
1815  replace_aseed(const std::vector<std::vector<M> >& aseed, casadi_int npar) const {
1816  std::vector<std::vector<M> > r(aseed.size());
1817  for (casadi_int d=0; d<r.size(); ++d) r[d] = replace_res(aseed[d], npar);
1818  return r;
1819  }
1820 
1821  template<typename M>
1822  std::map<std::string, M> FunctionInternal::
1823  convert_arg(const std::vector<M>& arg) const {
1824  casadi_assert(arg.size()==n_in_, "Incorrect number of inputs: Expected "
1825  + str(n_in_) + ", got " + str(arg.size()));
1826  std::map<std::string, M> ret;
1827  for (casadi_int i=0;i<n_in_;++i) {
1828  ret[name_in_[i]] = arg[i];
1829  }
1830  return ret;
1831  }
1832 
1833  template<typename M>
1834  std::vector<M> FunctionInternal::
1835  convert_arg(const std::map<std::string, M>& arg) const {
1836  // Get default inputs
1837  std::vector<M> arg_v(n_in_);
1838  for (casadi_int i=0; i<arg_v.size(); ++i) {
1839  arg_v[i] = get_default_in(i);
1840  }
1841 
1842  // Assign provided inputs
1843  for (auto&& e : arg) {
1844  arg_v.at(index_in(e.first)) = e.second;
1845  }
1846 
1847  return arg_v;
1848  }
1849 
1850  template<typename M>
1851  std::map<std::string, M> FunctionInternal::
1852  convert_res(const std::vector<M>& res) const {
1853  casadi_assert(res.size()==n_out_, "Incorrect number of outputs: Expected "
1854  + str(n_out_) + ", got " + str(res.size()));
1855  std::map<std::string, M> ret;
1856  for (casadi_int i=0;i<n_out_;++i) {
1857  ret[name_out_[i]] = res[i];
1858  }
1859  return ret;
1860  }
1861 
1862  template<typename M>
1863  std::vector<M> FunctionInternal::
1864  convert_res(const std::map<std::string, M>& res) const {
1865  // Get default inputs
1866  std::vector<M> res_v(n_out_);
1867  for (casadi_int i=0; i<res_v.size(); ++i) {
1868  res_v[i] = std::numeric_limits<double>::quiet_NaN();
1869  }
1870 
1871  // Assign provided inputs
1872  for (auto&& e : res) {
1873  M a = e.second;
1874  res_v.at(index_out(e.first)) = a;
1875  }
1876  return res_v;
1877  }
1878 
1879  template<typename MatType>
1880  MatType FunctionInternal::ensure_stacked(const MatType& v, const Sparsity& sp, casadi_int n) {
1881  // Check dimensions
1882  if (v.size1() == sp.size1() && v.size2() == n * sp.size2()) {
1883  // Ensure that sparsity is a horizontal multiple of original input, or has no entries
1884  if (v.nnz() != 0 && !v.sparsity().is_stacked(sp, n)) {
1885  return project(v, repmat(sp, 1, n));
1886  }
1887  } else {
1888  // Correct empty sparsity
1889  casadi_assert_dev(v.is_empty());
1890  return MatType(sp.size1(), sp.size2() * n);
1891  }
1892  // No correction needed
1893  return v;
1894  }
1895 
1896 } // namespace casadi
1897 
1899 
1900 #endif // CASADI_FUNCTION_INTERNAL_HPP
casadi::Function::AuxOut
std::map< std::string, std::vector< std::string > > AuxOut
Definition: function.hpp:404
casadi::GenericMatrix< Matrix< Scalar > >::zeros
static Matrix< Scalar > zeros(casadi_int nrow=1, casadi_int ncol=1)
Create a dense matrix or a matrix with specified sparsity with all entries zero.
Definition: generic_matrix.hpp:1262
casadi
The casadi namespace.
Definition: archiver.hpp:32
casadi::MXVector
std::vector< MX > MXVector
Definition: mx.hpp:1006
casadi::SX
Matrix< SXElem > SX
Definition: sx_fwd.hpp:32
casadi::Dict
GenericType::Dict Dict
C++ equivalent of Python's dict or MATLAB's struct.
Definition: generic_type.hpp:258
casadi::DM
Matrix< double > DM
Definition: dm_fwd.hpp:33
casadi::UnifiedReturnStatus
UnifiedReturnStatus
Definition: casadi_types.hpp:55