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  // Shared sequence number for input, output and instruction dumps.
86  casadi_int dump_id;
87  FunctionMemory() : stats_available(false), dump_id(-1) {}
88  };
89 
96  class CASADI_EXPORT ProtoFunction : public SharedObjectInternal {
97  public:
101  ProtoFunction(const std::string& name);
102 
106  ~ProtoFunction() override = 0;
107 
113  void construct(const Dict& opts);
114 
116 
119  static const Options options_;
120  virtual const Options& get_options() const { return options_;}
122 
124  virtual Dict generate_options(const std::string& target) const;
125 
129  void print_options(std::ostream &stream) const;
130 
134  void print_option(const std::string &name, std::ostream &stream) const;
135 
139  bool has_option(const std::string &option_name) const;
140 
144  virtual void change_option(const std::string& option_name, const GenericType& option_value);
145 
154  virtual void init(const Dict& opts);
155 
162  virtual void finalize();
163 
165  int checkout() const;
166 
168  void release(int mem) const;
169 
171  void* memory(int ind) const;
172 
174  bool has_memory(int ind) const;
175 
181  virtual void check_mem_count(casadi_int n) const { }
182 
186  virtual void* alloc_mem() const { return new ProtoFunctionMemory(); }
187 
191  virtual int init_mem(void* mem) const;
192 
196  virtual void free_mem(void *mem) const { delete static_cast<ProtoFunctionMemory*>(mem); }
197 
199  virtual Dict get_stats(void* mem) const;
200 
204  void clear_mem();
205 
209  void print(const char* fmt, ...) const;
210 
214  void sprint(char* buf, size_t buf_sz, const char* fmt, ...) const;
215 
219  void format_time(char* buffer, double time) const;
220 
224  void print_time(const std::map<std::string, FStats>& fstats) const;
225 
229  void serialize(SerializingStream &s) const;
230 
234  virtual void serialize_body(SerializingStream &s) const;
238  virtual void serialize_type(SerializingStream &s) const {}
239 
243  virtual std::string serialize_base_function() const {
244  return class_name();
245  }
246 
248  std::string name_;
249 
251  bool verbose_;
252 
253  // Print timing statistics
254  bool print_time_;
255 
256  // Print timing statistics
257  bool record_time_;
258 
260  bool regularity_check_;
261 
263  bool error_on_fail_;
264 
265  protected:
269  explicit ProtoFunction(DeserializingStream& s);
270 
271 #ifdef CASADI_WITH_THREAD
273  mutable std::mutex mtx_;
274 #endif // CASADI_WITH_THREAD
275 
276  private:
278  mutable std::vector<void*> mem_;
279 
281  mutable std::stack<int> unused_;
282  };
283 
290  class CASADI_EXPORT FunctionInternal : public ProtoFunction {
291  friend class Function;
292  public:
296  FunctionInternal(const std::string& name);
297 
301  ~FunctionInternal() override = 0;
302 
306  virtual std::string getAdaptorSolverName() const { return ""; }
307 
309 
312  static const Options options_;
313  const Options& get_options() const override { return options_;}
315 
317  Dict generate_options(const std::string& target) const override;
318 
322  void change_option(const std::string& option_name, const GenericType& option_value) override;
323 
327  void reset_dump_count();
328 
332  void init(const Dict& opts) override;
333 
337  void finalize() override;
338 
342  void* alloc_mem() const override { return new FunctionMemory(); }
343 
347  void free_mem(void *mem) const override { delete static_cast<FunctionMemory*>(mem); }
348 
350  Dict get_stats(void* mem) const override;
351 
355  Function self() const { return shared_from_this<Function>();}
356 
357  // Factory
358  virtual Function factory(const std::string& name,
359  const std::vector<std::string>& s_in,
360  const std::vector<std::string>& s_out,
361  const Function::AuxOut& aux,
362  const Dict& opts) const;
363 
364  // Get list of dependency functions
365  virtual std::vector<std::string> get_function() const;
366 
367  // Get a dependency function
368  virtual const Function& get_function(const std::string &name) const;
369 
370  // Check if a particular dependency exists
371  virtual bool has_function(const std::string& fname) const {return false;}
372 
373  // Add embedded function to map, helper function
374  void add_embedded(std::map<FunctionInternal*, std::pair<Function, size_t> >& all_fun,
375  const Function& dep, casadi_int max_depth) const;
376 
377  // Get all embedded functions, recursively
378  virtual void find(std::map<FunctionInternal*, std::pair<Function, size_t> >& all_fun,
379  casadi_int max_depth) const;
380 
389  virtual std::vector<bool> which_depends(const std::string& s_in,
390  const std::vector<std::string>& s_out,
391  casadi_int order, bool tr=false) const;
392 
399  virtual Function simplify_passes(
400  const std::vector<std::pair<std::string, casadi_int> >& tasks) const;
401 
403 
406  virtual bool has_spfwd() const { return false;}
407  virtual bool has_sprev() const { return false;}
409 
411 
414  int eval_gen(const double** arg, double** res, casadi_int* iw, double* w, void* mem,
415  bool always_inline, bool never_inline) const;
416  virtual int eval(const double** arg, double** res, casadi_int* iw, double* w, void* mem) const;
417 
419 
423  virtual int eval_sx(const SXElem** arg, SXElem** res,
424  casadi_int* iw, SXElem* w, void* mem, bool always_inline, bool never_inline) const;
425 
429  virtual void eval_mx(const MXVector& arg, MXVector& res,
430  bool always_inline, bool never_inline) const;
431 
433 
436  virtual std::vector<DM> eval_dm(const std::vector<DM>& arg) const;
437  virtual bool has_eval_dm() const { return false;}
439 
441 
444  int eval_gen(const SXElem** arg, SXElem** res, casadi_int* iw, SXElem* w, void* mem,
445  bool always_inline, bool never_inline) const {
446  return eval_sx(arg, res, iw, w, mem, always_inline, never_inline);
447  }
448  int eval_gen(const bvec_t** arg, bvec_t** res, casadi_int* iw, bvec_t* w, void* mem,
449  bool always_inline, bool never_inline) const {
450  return sp_forward(arg, res, iw, w, mem);
451  }
453 
455 
458  void call_gen(const MXVector& arg, MXVector& res, casadi_int npar,
459  bool always_inline, bool never_inline) const;
460 
461  template<typename D>
462  void call_gen(const std::vector<Matrix<D> >& arg, std::vector<Matrix<D> >& res,
463  casadi_int npar, bool always_inline, bool never_inline) const;
465 
469  template<typename M>
470  void call(const std::vector<M>& arg, std::vector<M>& res,
471  bool always_inline, bool never_inline) const;
472 
474 
479  static bool check_mat(const Sparsity& arg, const Sparsity& inp, casadi_int& npar);
481 
483 
491  template<typename M>
492  void check_arg(const std::vector<M>& arg, casadi_int& npar) const;
494 
496 
504  template<typename M>
505  void check_res(const std::vector<M>& res, casadi_int& npar) const;
507 
516  template<typename M> bool
517  matching_arg(const std::vector<M>& arg, casadi_int& npar) const;
518 
527  template<typename M> bool
528  matching_res(const std::vector<M>& arg, casadi_int& npar) const;
529 
533  template<typename M> std::vector<M>
534  replace_arg(const std::vector<M>& arg, casadi_int npar) const;
535 
539  template<typename M> std::vector<M>
540  project_arg(const std::vector<M>& arg, casadi_int npar) const;
541 
545  template<typename M> std::vector<M>
546  project_res(const std::vector<M>& arg, casadi_int npar) const;
547 
551  template<typename M> std::vector<M>
552  replace_res(const std::vector<M>& res, casadi_int npar) const;
553 
557  template<typename M> std::vector<std::vector<M>>
558  replace_fseed(const std::vector<std::vector<M>>& fseed, casadi_int npar) const;
559 
563  template<typename M> std::vector<std::vector<M>>
564  replace_aseed(const std::vector<std::vector<M>>& aseed, casadi_int npar) const;
565 
570  template<typename M>
571  std::map<std::string, M> convert_arg(const std::vector<M>& arg) const;
572  template<typename M>
573  std::vector<M> convert_arg(const std::map<std::string, M>& arg) const;
574  template<typename M>
575  std::map<std::string, M> convert_res(const std::vector<M>& res) const;
576  template<typename M>
577  std::vector<M> convert_res(const std::map<std::string, M>& res) const;
579 
584  std::vector<double> nz_in(const std::vector<DM>& arg) const;
585  std::vector<double> nz_out(const std::vector<DM>& res) const;
586  std::vector<DM> nz_in(const std::vector<double>& arg) const;
587  std::vector<DM> nz_out(const std::vector<double>& res) const;
589 
591 
594  virtual void call_forward(const std::vector<MX>& arg, const std::vector<MX>& res,
595  const std::vector<std::vector<MX> >& fseed,
596  std::vector<std::vector<MX> >& fsens,
597  bool always_inline, bool never_inline) const;
598  virtual void call_forward(const std::vector<SX>& arg, const std::vector<SX>& res,
599  const std::vector<std::vector<SX> >& fseed,
600  std::vector<std::vector<SX> >& fsens,
601  bool always_inline, bool never_inline) const;
603 
605 
608  virtual void call_reverse(const std::vector<MX>& arg, const std::vector<MX>& res,
609  const std::vector<std::vector<MX> >& aseed,
610  std::vector<std::vector<MX> >& asens,
611  bool always_inline, bool never_inline) const;
612  virtual void call_reverse(const std::vector<SX>& arg, const std::vector<SX>& res,
613  const std::vector<std::vector<SX> >& aseed,
614  std::vector<std::vector<SX> >& asens,
615  bool always_inline, bool never_inline) const;
617 
621  std::vector<MX> mapsum_mx(const std::vector<MX > &arg, const std::string& parallelization);
622 
626  virtual bool uses_output() const {return false;}
627 
629 
632  Function jacobian() const;
633  virtual bool has_jacobian() const { return false;}
634  virtual Function get_jacobian(const std::string& name,
635  const std::vector<std::string>& inames,
636  const std::vector<std::string>& onames,
637  const Dict& opts) const;
639 
641 
645  Sparsity& jac_sparsity(casadi_int oind, casadi_int iind, bool compact, bool symmetric) const;
646  virtual bool has_jac_sparsity(casadi_int oind, casadi_int iind) const { return false;}
647  virtual Sparsity get_jac_sparsity(casadi_int oind, casadi_int iind, bool symmetric) const;
649 
651  static std::string forward_name(const std::string& fcn, casadi_int nfwd) {
652  return "fwd" + str(nfwd) + "_" + fcn;
653  }
654 
656  std::string diff_prefix(const std::string& prefix) const;
657 
659 
666  Function forward(casadi_int nfwd) const;
667  virtual bool has_forward(casadi_int nfwd) const { return false;}
668  virtual Function get_forward(casadi_int nfwd, const std::string& name,
669  const std::vector<std::string>& inames,
670  const std::vector<std::string>& onames,
671  const Dict& opts) const;
673 
675  static std::string reverse_name(const std::string& fcn, casadi_int nadj) {
676  return "adj" + str(nadj) + "_" + fcn;
677  }
678 
680 
687  Function reverse(casadi_int nadj) const;
688  virtual bool has_reverse(casadi_int nadj) const { return false;}
689  virtual Function get_reverse(casadi_int nadj, const std::string& name,
690  const std::vector<std::string>& inames,
691  const std::vector<std::string>& onames,
692  const Dict& opts) const;
694 
698  template<typename MatType>
699  static MatType ensure_stacked(const MatType& v, const Sparsity& sp, casadi_int n);
700 
704  virtual Function slice(const std::string& name, const std::vector<casadi_int>& order_in,
705  const std::vector<casadi_int>& order_out, const Dict& opts) const;
706 
710  virtual const Function& oracle() const;
711 
715  bool has_derivative() const;
716 
720  virtual double ad_weight() const;
721 
727  virtual double sp_weight() const;
728 
730 
733  virtual const SX sx_in(casadi_int ind) const;
734  virtual const SX sx_out(casadi_int ind) const;
735  virtual const std::vector<SX> sx_in() const;
736  virtual const std::vector<SX> sx_out() const;
737  virtual const MX mx_in(casadi_int ind) const;
738  virtual const MX mx_out(casadi_int ind) const;
739  virtual const std::vector<MX> mx_in() const;
740  virtual const std::vector<MX> mx_out() const;
741  const DM dm_in(casadi_int ind) const;
742  const DM dm_out(casadi_int ind) const;
743  const std::vector<DM> dm_in() const;
744  const std::vector<DM> dm_out() const;
746 
748  virtual std::vector<MX> free_mx() const;
749 
751  virtual std::vector<SX> free_sx() const;
752 
756  virtual bool has_free() const { return false;}
757 
761  virtual void generate_lifted(Function& vdef_fcn, Function& vinit_fcn) const;
762 
766  virtual casadi_int n_instructions() const;
767 
771  virtual casadi_int instruction_id(casadi_int k) const;
772 
776  virtual std::vector<casadi_int> instruction_input(casadi_int k) const;
777 
781  virtual double instruction_constant(casadi_int k) const;
782 
786  virtual std::vector<casadi_int> instruction_output(casadi_int k) const;
787 
791  virtual casadi_int n_nodes() const;
792 
796  virtual MX instruction_MX(casadi_int k) const;
797 
801  virtual SX instructions_sx() const;
802 
807  Function wrap(const std::string& name) const;
808  Function wrap() const;
810 
815  Function wrap_as_needed(const std::string& name, const Dict& opts) const;
816  Function wrap_as_needed(const Dict& opts) const;
818 
822  Dict cache() const;
823 
827  bool incache(const std::string& fname, Function& f, const std::string& suffix="") const;
828 
832  void tocache(const Function& f, const std::string& suffix="") const;
833 
834 
838  void tocache_if_missing(Function& f, const std::string& suffix="") const;
839 
843  void codegen(CodeGenerator& g, const std::string& fname) const;
844 
848  void codegen_meta(CodeGenerator& g) const;
849 
853  void codegen_sparsities(CodeGenerator& g) const;
854 
858  virtual std::string codegen_name(const CodeGenerator& g, bool ns=true) const;
859 
863  std::string codegen_mem(CodeGenerator& g, const std::string& index="mem") const;
864 
868  virtual void codegen_incref(CodeGenerator& g) const;
869 
873  virtual void codegen_decref(CodeGenerator& g) const;
874 
878  virtual void codegen_alloc_mem(CodeGenerator& g) const;
879 
883  virtual void codegen_init_mem(CodeGenerator& g) const;
884 
888  virtual void codegen_free_mem(CodeGenerator& g) const {}
889 
893  virtual void codegen_checkout(CodeGenerator& g) const;
894 
898  virtual void codegen_release(CodeGenerator& g) const;
899 
903  std::string signature(const std::string& fname) const;
904 
908  std::string signature_unrolled(const std::string& fname) const;
909 
913  virtual void codegen_declarations(CodeGenerator& g) const;
914 
918  virtual void codegen_body(CodeGenerator& g) const;
919 
923  virtual std::string codegen_mem_type() const { return ""; }
924 
928  virtual bool codegen_needs_mem() const { return false; }
929 
935  virtual bool codegen_mem_is_opaque() const { return false; }
936 
940  virtual std::string generate_dependencies(const std::string& fname, const Dict& opts) const;
941 
945  virtual bool has_codegen() const { return false;}
946 
950  virtual void jit_dependencies(const std::string& fname) {}
951 
955  static std::string get_jit_directory(const Dict& jit_options);
956 
960  virtual void export_code(const std::string& lang,
961  std::ostream &stream, const Dict& options) const;
962 
966  void serialize_type(SerializingStream &s) const override;
967 
971  void serialize_body(SerializingStream &s) const override;
972 
976  void disp(std::ostream& stream, bool more) const override;
977 
981  virtual void disp_more(std::ostream& stream) const {}
982 
986  std::string definition() const;
987 
991  void print_dimensions(std::ostream &stream) const;
992 
996  virtual std::vector<std::string> get_free() const;
997 
1001  void get_partition(casadi_int iind, casadi_int oind, Sparsity& D1, Sparsity& D2,
1002  bool compact, bool symmetric,
1003  bool allow_forward, bool allow_reverse) const;
1004 
1006 
1009  casadi_int nnz_in() const;
1010  casadi_int nnz_in(casadi_int ind) const { return sparsity_in(ind).nnz(); }
1011  casadi_int nnz_out() const;
1012  casadi_int nnz_out(casadi_int ind) const { return sparsity_out(ind).nnz(); }
1014 
1016 
1019  casadi_int numel_in() const;
1020  casadi_int numel_in(casadi_int ind) const { return sparsity_in(ind).numel(); }
1021  casadi_int numel_out(casadi_int ind) const { return sparsity_out(ind).numel(); }
1022  casadi_int numel_out() const;
1024 
1026 
1029  casadi_int size1_in(casadi_int ind) const { return sparsity_in(ind).size1(); }
1030  casadi_int size2_in(casadi_int ind) const { return sparsity_in(ind).size2(); }
1031  casadi_int size1_out(casadi_int ind) const { return sparsity_out(ind).size1(); }
1032  casadi_int size2_out(casadi_int ind) const { return sparsity_out(ind).size2(); }
1033  std::pair<casadi_int, casadi_int> size_in(casadi_int ind) const {
1034  return sparsity_in(ind).size();
1035  }
1036  std::pair<casadi_int, casadi_int> size_out(casadi_int ind) const {
1037  return sparsity_out(ind).size();
1038  }
1040 
1042 
1045  const Sparsity& sparsity_in(casadi_int ind) const { return sparsity_in_.at(ind); }
1046  const Sparsity& sparsity_out(casadi_int ind) const { return sparsity_out_.at(ind); }
1048 
1050 
1053  bool all_scalar() const;
1054 
1056  virtual bool jac_is_symm(casadi_int oind, casadi_int iind) const;
1057 
1059  Sparsity to_compact(casadi_int oind, casadi_int iind, const Sparsity& sp) const;
1060 
1062  Sparsity from_compact(casadi_int oind, casadi_int iind, const Sparsity& sp) const;
1063 
1065  template<bool fwd>
1066  Sparsity get_jac_sparsity_gen(casadi_int oind, casadi_int iind) const;
1067 
1069  Sparsity get_jac_sparsity_hierarchical(casadi_int oind, casadi_int iind) const;
1070 
1074  Sparsity get_jac_sparsity_hierarchical_symm(casadi_int oind, casadi_int iind) const;
1075 
1077  virtual std::vector<MX> symbolic_output(const std::vector<MX>& arg) const;
1078 
1080 
1083  virtual size_t get_n_in();
1084  virtual size_t get_n_out();
1086 
1087 
1089 
1092  virtual std::string get_name_in(casadi_int i);
1093  virtual std::string get_name_out(casadi_int i);
1095 
1099  virtual double get_default_in(casadi_int ind) const {
1100  return 0;
1101  }
1102 
1106  virtual double get_max_in(casadi_int ind) const {
1107  return inf;
1108  }
1109 
1113  virtual double get_min_in(casadi_int ind) const {
1114  return -inf;
1115  }
1116 
1117  virtual std::vector<double> get_nominal_in(casadi_int ind) const {
1118  return std::vector<double>(nnz_in(ind), 1.);
1119  }
1120 
1121  virtual std::vector<double> get_nominal_out(casadi_int ind) const {
1122  return std::vector<double>(nnz_out(ind), 1.);
1123  }
1124 
1128  virtual double get_reltol() const {
1129  return eps;
1130  }
1131 
1135  virtual double get_abstol() const {
1136  return eps;
1137  }
1138 
1142  virtual Sparsity get_sparsity_in(casadi_int i);
1143 
1147  virtual Sparsity get_sparsity_out(casadi_int i);
1148 
1152  virtual bool get_diff_in(casadi_int i) { return true; }
1153 
1157  virtual bool get_diff_out(casadi_int i) { return true; }
1158 
1162  casadi_int index_in(const std::string &name) const {
1163  for (casadi_int i=0; i<name_in_.size(); ++i) {
1164  if (name_in_[i]==name) return i;
1165  }
1166  casadi_error("FunctionInternal::index_in: could not find entry \""
1167  + name + "\". Available names are: " + str(name_in_) + ".");
1168  return -1;
1169  }
1170 
1174  casadi_int index_out(const std::string &name) const {
1175  for (casadi_int i=0; i<name_out_.size(); ++i) {
1176  if (name_out_[i]==name) return i;
1177  }
1178  casadi_error("FunctionInternal::index_out: could not find entry \""
1179  + name + "\". Available names are: " + str(name_out_) + ".");
1180  return -1;
1181  }
1182 
1186  virtual int sp_forward(const bvec_t** arg, bvec_t** res,
1187  casadi_int* iw, bvec_t* w, void* mem) const;
1188 
1196  virtual int eval_activity(const bvec_t** arg, bvec_t** res,
1197  casadi_int* iw, bvec_t* w, void* mem) const;
1198 
1202  virtual int sp_forward_block(const bvec_t** arg, bvec_t** res,
1203  casadi_int* iw, bvec_t* w, void* mem, casadi_int oind, casadi_int iind) const;
1204 
1208  virtual int sp_reverse(bvec_t** arg, bvec_t** res, casadi_int* iw, bvec_t* w, void* mem) const;
1209 
1213  void sz_work(size_t& sz_arg, size_t& sz_res, size_t& sz_iw, size_t& sz_w) const;
1214 
1218  size_t sz_arg() const { return sz_arg_per_ + sz_arg_tmp_;}
1219 
1223  size_t sz_res() const { return sz_res_per_ + sz_res_tmp_;}
1224 
1228  size_t sz_iw() const { return sz_iw_per_ + sz_iw_tmp_;}
1229 
1233  size_t sz_w() const { return sz_w_per_ + sz_w_tmp_;}
1234 
1239  virtual size_t codegen_sz_arg(const CodeGenerator& g) const;
1240  virtual size_t codegen_sz_res(const CodeGenerator& g) const;
1241  virtual size_t codegen_sz_iw(const CodeGenerator& g) const;
1242  virtual size_t codegen_sz_w(const CodeGenerator& g) const;
1244 
1248  void alloc_arg(size_t sz_arg, bool persistent=false);
1249 
1253  void alloc_res(size_t sz_res, bool persistent=false);
1254 
1258  void alloc_iw(size_t sz_iw, bool persistent=false);
1259 
1263  void alloc_w(size_t sz_w, bool persistent=false);
1264 
1268  void alloc(const Function& f, bool persistent=false, int num_threads=1);
1269 
1273  virtual void set_work(void* mem, const double**& arg, double**& res,
1274  casadi_int*& iw, double*& w) const {}
1275 
1279  virtual void set_temp(void* mem, const double** arg, double** res,
1280  casadi_int* iw, double* w) const {}
1281 
1285  void setup(void* mem, const double** arg, double** res, casadi_int* iw, double* w) const;
1286 
1288 
1291  virtual bool fwdViaJac(casadi_int nfwd) const;
1292  virtual bool adjViaJac(casadi_int nadj) const;
1294 
1296  virtual Dict info() const;
1297 
1301  Function map(casadi_int n, const std::string& parallelization) const;
1302 
1306  void generate_in(const std::string& fname, const double** arg) const;
1307  void generate_out(const std::string& fname, double** res) const;
1308 
1309  bool always_inline_, never_inline_;
1310 
1312  size_t n_in_, n_out_;
1313 
1315  std::vector<bool> is_diff_in_, is_diff_out_;
1316 
1318  std::vector<Sparsity> sparsity_in_, sparsity_out_;
1319 
1321  std::vector<std::string> name_in_, name_out_;
1322 
1326  bool jit_;
1327 
1331  bool jit_cleanup_;
1332 
1336  std::string jit_serialize_;
1337 
1341  std::string jit_name_;
1342  std::string jit_directory_;
1343  std::string jit_base_name_;
1344 
1348  bool jit_temp_suffix_;
1349 
1353  eval_t eval_;
1354 
1358  casadi_checkout_t checkout_;
1359 
1363  casadi_release_t release_;
1364 
1368  signal_t incref_, decref_;
1369 
1373  Dict stats_;
1374 
1378  bool has_refcount_;
1379 
1383  bool has_refcount_in_deps_;
1384 
1388  Dict cache_init_;
1389 
1391  mutable WeakCache<std::string, Function> cache_;
1392 
1394  mutable std::vector<Sparsity> jac_sparsity_[2];
1395 
1396 #ifdef CASADI_WITH_THREADSAFE_SYMBOLICS
1398  mutable std::mutex jac_sparsity_mtx_;
1399 #endif // CASADI_WITH_THREADSAFE_SYMBOLICS
1400 
1402  Function derivative_of_;
1403 
1405  void* user_data_;
1406 
1408  std::string compiler_plugin_;
1409  Importer compiler_;
1410  Dict jit_options_;
1411 
1413  double jac_penalty_;
1414 
1415  // Types of derivative calculation permitted
1416  bool enable_forward_, enable_reverse_, enable_jacobian_, enable_fd_;
1417  bool enable_forward_op_, enable_reverse_op_, enable_jacobian_op_, enable_fd_op_;
1418 
1420  double ad_weight_, ad_weight_sp_;
1421 
1423  casadi_int max_num_dir_;
1424 
1426  bool inputs_check_;
1427 
1428  // Finite difference step
1429  Dict fd_options_;
1430 
1431  // Finite difference step size
1432  double fd_step_;
1433 
1434  // Finite difference method
1435  std::string fd_method_;
1436 
1437  // Print input/output
1438  bool print_in_;
1439  bool print_out_;
1440 
1441  // Print canonical form
1442  bool print_canonical_;
1443 
1444  // Warn when number of inputs or outputs exceed this value
1445  casadi_int max_io_;
1446 
1447  // Dump input/output
1448  bool dump_in_, dump_out_, dump_;
1449 
1450  // Directory to dump to
1451  std::string dump_dir_;
1452 
1453  // Format to dump with
1454  std::string dump_format_;
1455 
1456  // Forward/reverse/Jacobian options
1457  Dict forward_options_, reverse_options_, jacobian_options_, der_options_;
1458 
1459  // Store a reference to a custom Jacobian
1460  Function custom_jacobian_;
1461 
1462  // Registered functions
1463  std::vector<Function> registered_functions_;
1464 
1465  // Counter for unique names for dumping inputs and output
1466 #ifdef CASADI_WITH_THREAD
1467  mutable std::atomic<casadi_int> dump_count_;
1468 #else
1469  mutable casadi_int dump_count_;
1470 #endif // CASADI_WITH_THREAD
1471 
1475  virtual bool is_a(const std::string& type, bool recursive) const;
1476 
1480  virtual void merge(const std::vector<MX>& arg,
1481  std::vector<MX>& subs_from, std::vector<MX>& subs_to) const;
1482 
1486  template<typename MatType>
1487  static bool purgable(const std::vector<MatType>& seed);
1488 
1492  template<typename MatType>
1493  std::vector<std::vector<MatType> >
1494  fwd_seed(casadi_int nfwd) const;
1495 
1499  template<typename MatType>
1500  std::vector<std::vector<MatType> >
1501  symbolicAdjSeed(casadi_int nadj, const std::vector<MatType>& v) const;
1502 
1503  static std::string string_from_UnifiedReturnStatus(UnifiedReturnStatus status);
1504 
1508  explicit FunctionInternal(DeserializingStream& e);
1509 
1513  static Function deserialize(DeserializingStream& s);
1514  static std::map<std::string, ProtoFunction* (*)(DeserializingStream&)> deserialize_map;
1515 
1519  void print_in(std::ostream &stream, const double** arg, bool truncate) const;
1520 
1524  void print_out(std::ostream &stream, double** res, bool truncate) const;
1525 
1529  static void print_canonical(std::ostream &stream, const Sparsity& sp, const double* nz);
1530 
1534  static void print_canonical(std::ostream &stream, casadi_int sz, const double* nz);
1535 
1539  static void print_canonical(std::ostream &stream, double a);
1540 
1541  protected:
1545  void set_jac_sparsity(casadi_int oind, casadi_int iind, const Sparsity& sp);
1546 
1547  std::unique_ptr<std::ostream> open_trace(const double** arg, casadi_int dump_id) const;
1548  void finish_trace(std::ostream& trace, double** res, int ret) const;
1549  static void trace_values(std::ostream& trace, const double* values, casadi_int nnz);
1550 
1551  private:
1552  // @{
1554  casadi_int get_dump_id() const;
1555  void dump_in(casadi_int id, const double** arg) const;
1556  void dump_out(casadi_int id, double** res) const;
1557  void dump() const;
1558  // @}
1559 
1563  size_t sz_arg_per_, sz_res_per_, sz_iw_per_, sz_w_per_;
1564 
1568  size_t sz_arg_tmp_, sz_res_tmp_, sz_iw_tmp_, sz_w_tmp_;
1569  };
1570 
1571  // Template implementations
1572  template<typename MatType>
1573  bool FunctionInternal::purgable(const std::vector<MatType>& v) {
1574  for (auto i=v.begin(); i!=v.end(); ++i) {
1575  if (!i->is_zero()) return false;
1576  }
1577  return true;
1578  }
1579 
1580  template<typename MatType>
1581  std::vector<std::vector<MatType> >
1582  FunctionInternal::
1583  fwd_seed(casadi_int nfwd) const {
1584  std::vector<std::vector<MatType>> fseed(nfwd);
1585  for (casadi_int dir=0; dir<nfwd; ++dir) {
1586  fseed[dir].resize(n_in_);
1587  for (casadi_int iind=0; iind<n_in_; ++iind) {
1588  std::string n = "f" + str(dir) + "_" + name_in_[iind];
1589  Sparsity sp = is_diff_in_[iind] ? sparsity_in(iind) : Sparsity(size_in(iind));
1590  fseed[dir][iind] = MatType::sym(n, sp);
1591  }
1592  }
1593  return fseed;
1594  }
1595 
1596  template<typename MatType>
1597  std::vector<std::vector<MatType> >
1598  FunctionInternal::
1599  symbolicAdjSeed(casadi_int nadj, const std::vector<MatType>& v) const {
1600  std::vector<std::vector<MatType> > aseed(nadj, v);
1601  for (casadi_int dir=0; dir<nadj; ++dir) {
1602  // Replace symbolic inputs
1603  casadi_int oind=0;
1604  for (typename std::vector<MatType>::iterator i=aseed[dir].begin();
1605  i!=aseed[dir].end();
1606  ++i, ++oind) {
1607  // Name of the adjoint seed
1608  std::stringstream ss;
1609  ss << "a";
1610  if (nadj>1) ss << dir << "_";
1611  ss << oind;
1612 
1613  // Save to matrix
1614  *i = MatType::sym(ss.str(), is_diff_out_[oind] ? i->sparsity() : Sparsity(i->size()));
1615 
1616  }
1617  }
1618  return aseed;
1619  }
1620 
1621  template<typename M>
1622  void FunctionInternal::call(const std::vector<M>& arg, std::vector<M>& res,
1623  bool always_inline, bool never_inline) const {
1624  // If all inputs are scalar ...
1625  if (all_scalar()) {
1626  // ... and some arguments are matrix-valued with matching dimensions ...
1627  bool matrix_call = false;
1628  std::pair<casadi_int, casadi_int> sz;
1629  for (auto&& a : arg) {
1630  if (!a.is_scalar() && !a.is_empty()) {
1631  if (!matrix_call) {
1632  // Matrix call
1633  matrix_call = true;
1634  sz = a.size();
1635  } else if (a.size()!=sz) {
1636  // Not same dimensions
1637  matrix_call = false;
1638  break;
1639  }
1640  }
1641  }
1642 
1643  // ... then, call multiple times
1644  if (matrix_call) {
1645  // Start with zeros
1646  res.resize(n_out_);
1647  M z = M::zeros(sz);
1648  for (auto&& a : res) a = z;
1649  // Call multiple times
1650  std::vector<M> arg1 = arg, res1;
1651  for (casadi_int c=0; c<sz.second; ++c) {
1652  for (casadi_int r=0; r<sz.first; ++r) {
1653  // Get scalar arguments
1654  for (casadi_int i=0; i<arg.size(); ++i) {
1655  if (arg[i].size()==sz) arg1[i] = arg[i](r, c);
1656  }
1657  // Call recursively with scalar arguments
1658  call(arg1, res1, always_inline, never_inline);
1659  // Get results
1660  casadi_assert_dev(res.size() == res1.size());
1661  for (casadi_int i=0; i<res.size(); ++i) res[i](r, c) = res1[i];
1662  }
1663  }
1664  // All elements assigned
1665  return;
1666  }
1667  }
1668 
1669  // Check if inputs need to be replaced
1670  casadi_int npar = 1;
1671  if (!matching_arg(arg, npar)) {
1672  return call(replace_arg(arg, npar), res, always_inline, never_inline);
1673  }
1674 
1675  // Call the type-specific method
1676  call_gen(arg, res, npar, always_inline, never_inline);
1677  }
1678 
1679  template<typename M>
1680  std::vector<M> FunctionInternal::
1681  project_arg(const std::vector<M>& arg, casadi_int npar) const {
1682  casadi_assert_dev(arg.size()==n_in_);
1683 
1684  // Which arguments require mapped evaluation
1685  std::vector<bool> mapped(n_in_);
1686  for (casadi_int i=0; i<n_in_; ++i) {
1687  mapped[i] = arg[i].size2()!=size2_in(i);
1688  }
1689 
1690  // Check if matching input sparsity
1691  std::vector<bool> matching(n_in_);
1692  bool any_mismatch = false;
1693  for (casadi_int i=0; i<n_in_; ++i) {
1694  if (mapped[i]) {
1695  matching[i] = arg[i].sparsity().is_stacked(sparsity_in(i), npar);
1696  } else {
1697  matching[i] = arg[i].sparsity()==sparsity_in(i);
1698  }
1699  any_mismatch = any_mismatch || !matching[i];
1700  }
1701 
1702  // Correct input sparsity
1703  if (any_mismatch) {
1704  std::vector<M> arg2(arg);
1705  for (casadi_int i=0; i<n_in_; ++i) {
1706  if (!matching[i]) {
1707  if (mapped[i]) {
1708  arg2[i] = project(arg2[i], repmat(sparsity_in(i), 1, npar));
1709  } else {
1710  arg2[i] = project(arg2[i], sparsity_in(i));
1711  }
1712  }
1713  }
1714  return arg2;
1715  }
1716  return arg;
1717  }
1718 
1719  template<typename M>
1720  std::vector<M> FunctionInternal::
1721  project_res(const std::vector<M>& arg, casadi_int npar) const {
1722  return arg;
1723  }
1724 
1725  template<typename D>
1726  void FunctionInternal::
1727  call_gen(const std::vector<Matrix<D> >& arg, std::vector<Matrix<D> >& res,
1728  casadi_int npar, bool always_inline, bool never_inline) const {
1729  std::vector< Matrix<D> > arg2 = project_arg(arg, npar);
1730 
1731  // Which arguments require mapped evaluation
1732  std::vector<bool> mapped(n_in_);
1733  for (casadi_int i=0; i<n_in_; ++i) {
1734  mapped[i] = arg[i].size2()!=size2_in(i);
1735  }
1736 
1737  // Allocate results
1738  res.resize(n_out_);
1739  for (casadi_int i=0; i<n_out_; ++i) {
1740  if (!res[i].sparsity().is_stacked(sparsity_out(i), npar)) {
1741  res[i] = Matrix<D>::zeros(repmat(sparsity_out(i), 1, npar));
1742  }
1743  }
1744 
1745  // Allocate temporary memory if needed
1746  std::vector<casadi_int> iw_tmp(sz_iw());
1747  std::vector<D> w_tmp(sz_w());
1748 
1749  // Get pointers to input arguments
1750  std::vector<const D*> argp(sz_arg());
1751  for (casadi_int i=0; i<n_in_; ++i) argp[i]=get_ptr(arg2[i]);
1752 
1753  // Get pointers to output arguments
1754  std::vector<D*> resp(sz_res());
1755  for (casadi_int i=0; i<n_out_; ++i) resp[i]=get_ptr(res[i]);
1756 
1757  // For all parallel calls
1758  for (casadi_int p=0; p<npar; ++p) {
1759  // Call memory-less
1760  if (eval_gen(get_ptr(argp), get_ptr(resp),
1761  get_ptr(iw_tmp), get_ptr(w_tmp), memory(0),
1762  always_inline, never_inline)) {
1763  if (error_on_fail_) casadi_error("Evaluation failed");
1764  }
1765  // Update offsets
1766  if (p==npar-1) break;
1767  for (casadi_int i=0; i<n_in_; ++i) if (mapped[i]) argp[i] += nnz_in(i);
1768  for (casadi_int i=0; i<n_out_; ++i) resp[i] += nnz_out(i);
1769  }
1770  }
1771 
1772  template<typename M>
1773  void FunctionInternal::check_arg(const std::vector<M>& arg, casadi_int& npar) const {
1774  casadi_assert(arg.size()==n_in_, "Incorrect number of inputs: Expected "
1775  + str(n_in_) + ", got " + str(arg.size()));
1776  for (casadi_int i=0; i<n_in_; ++i) {
1777  if (!check_mat(arg[i].sparsity(), sparsity_in(i), npar)) {
1778  // Dimensions
1779  std::string d_arg = str(arg[i].size1()) + "-by-" + str(arg[i].size2());
1780  std::string d_in = str(size1_in(i)) + "-by-" + str(size2_in(i));
1781  std::string e = "Input " + str(i) + " (" + name_in_[i] + ") has mismatching shape. "
1782  "Got " + d_arg + ". Allowed dimensions, in general, are:\n"
1783  " - The input dimension N-by-M (here " + d_in + ")\n"
1784  " - A scalar, i.e. 1-by-1\n"
1785  " - M-by-N if N=1 or M=1 (i.e. a transposed vector)\n"
1786  " - N-by-M1 if K*M1=M for some K (argument repeated horizontally)\n";
1787  if (npar!=-1) {
1788  e += " - N-by-P*M, indicating evaluation with multiple arguments (P must be a "
1789  "multiple of " + str(npar) + " for consistency with previous inputs)";
1790  }
1791  casadi_error(e);
1792  }
1793  }
1794  }
1795 
1796  template<typename M>
1797  void FunctionInternal::check_res(const std::vector<M>& res, casadi_int& npar) const {
1798  casadi_assert(res.size()==n_out_, "Incorrect number of outputs: Expected "
1799  + str(n_out_) + ", got " + str(res.size()));
1800  for (casadi_int i=0; i<n_out_; ++i) {
1801  casadi_assert(check_mat(res[i].sparsity(), sparsity_out(i), npar),
1802  "Output " + str(i) + " (" + name_out_[i] + ") has mismatching shape. "
1803  "Expected " + str(size_out(i)) + ", got " + str(res[i].size()));
1804  }
1805  }
1806 
1807  template<typename M>
1808  bool FunctionInternal::matching_arg(const std::vector<M>& arg, casadi_int& npar) const {
1809  check_arg(arg, npar);
1810  for (casadi_int i=0; i<n_in_; ++i) {
1811  if (arg.at(i).size1()!=size1_in(i)) return false;
1812  if (arg.at(i).size2()!=size2_in(i) && arg.at(i).size2()!=npar*size2_in(i)) return false;
1813  }
1814  return true;
1815  }
1816 
1817  template<typename M>
1818  bool FunctionInternal::matching_res(const std::vector<M>& res, casadi_int& npar) const {
1819  check_res(res, npar);
1820  for (casadi_int i=0; i<n_out_; ++i) {
1821  if (res.at(i).size1()!=size1_out(i)) return false;
1822  if (res.at(i).size2()!=size2_out(i) && res.at(i).size2()!=npar*size2_out(i)) return false;
1823  }
1824  return true;
1825  }
1826 
1827  template<typename M>
1828  M replace_mat(const M& arg, const Sparsity& inp, casadi_int npar) {
1829  if (arg.size()==inp.size()) {
1830  // Matching dimensions already
1831  return arg;
1832  } else if (arg.is_scalar()) {
1833  // Scalar assign means set all
1834  return M(inp, arg);
1835  } else if (arg.is_vector() && inp.size()==std::make_pair(arg.size2(), arg.size1())) {
1836  // Transpose vector
1837  return arg.T();
1838  } else if (arg.size1()==inp.size1() && arg.size2()>0 && inp.size2()>0
1839  && inp.size2()%arg.size2()==0) {
1840  // Horizontal repmat
1841  return repmat(arg, 1, inp.size2()/arg.size2());
1842  } else if (npar!=-1 && arg.size1()==inp.size1() && arg.size2()>0 && inp.size2()>0
1843  && (npar*inp.size2())%arg.size2()==0) {
1844  // Multiple evaluation: grow argument horizontally to npar*inp columns
1845  return repmat(arg, 1, (npar*inp.size2())/arg.size2());
1846  } else {
1847  // Empty matrix means set zero (kept last so a 0-by-N argument is first given the
1848  // chance to be recognised as a parallel/repmat call above)
1849  casadi_assert_dev(arg.is_empty());
1850  return M(inp.size());
1851  }
1852  }
1853 
1854  template<typename M>
1855  std::vector<M> FunctionInternal::
1856  replace_arg(const std::vector<M>& arg, casadi_int npar) const {
1857  std::vector<M> r(arg.size());
1858  for (casadi_int i=0; i<r.size(); ++i) r[i] = replace_mat(arg[i], sparsity_in(i), npar);
1859  return r;
1860  }
1861 
1862  template<typename M>
1863  std::vector<M> FunctionInternal::
1864  replace_res(const std::vector<M>& res, casadi_int npar) const {
1865  std::vector<M> r(res.size());
1866  for (casadi_int i=0; i<r.size(); ++i) r[i] = replace_mat(res[i], sparsity_out(i), npar);
1867  return r;
1868  }
1869 
1870  template<typename M>
1871  std::vector<std::vector<M> > FunctionInternal::
1872  replace_fseed(const std::vector<std::vector<M> >& fseed, casadi_int npar) const {
1873  std::vector<std::vector<M> > r(fseed.size());
1874  for (casadi_int d=0; d<r.size(); ++d) r[d] = replace_arg(fseed[d], npar);
1875  return r;
1876  }
1877 
1878  template<typename M>
1879  std::vector<std::vector<M> > FunctionInternal::
1880  replace_aseed(const std::vector<std::vector<M> >& aseed, casadi_int npar) const {
1881  std::vector<std::vector<M> > r(aseed.size());
1882  for (casadi_int d=0; d<r.size(); ++d) r[d] = replace_res(aseed[d], npar);
1883  return r;
1884  }
1885 
1886  template<typename M>
1887  std::map<std::string, M> FunctionInternal::
1888  convert_arg(const std::vector<M>& arg) const {
1889  casadi_assert(arg.size()==n_in_, "Incorrect number of inputs: Expected "
1890  + str(n_in_) + ", got " + str(arg.size()));
1891  std::map<std::string, M> ret;
1892  for (casadi_int i=0;i<n_in_;++i) {
1893  ret[name_in_[i]] = arg[i];
1894  }
1895  return ret;
1896  }
1897 
1898  template<typename M>
1899  std::vector<M> FunctionInternal::
1900  convert_arg(const std::map<std::string, M>& arg) const {
1901  // Get default inputs
1902  std::vector<M> arg_v(n_in_);
1903  for (casadi_int i=0; i<arg_v.size(); ++i) {
1904  arg_v[i] = get_default_in(i);
1905  }
1906 
1907  // Assign provided inputs
1908  for (auto&& e : arg) {
1909  arg_v.at(index_in(e.first)) = e.second;
1910  }
1911 
1912  return arg_v;
1913  }
1914 
1915  template<typename M>
1916  std::map<std::string, M> FunctionInternal::
1917  convert_res(const std::vector<M>& res) const {
1918  casadi_assert(res.size()==n_out_, "Incorrect number of outputs: Expected "
1919  + str(n_out_) + ", got " + str(res.size()));
1920  std::map<std::string, M> ret;
1921  for (casadi_int i=0;i<n_out_;++i) {
1922  ret[name_out_[i]] = res[i];
1923  }
1924  return ret;
1925  }
1926 
1927  template<typename M>
1928  std::vector<M> FunctionInternal::
1929  convert_res(const std::map<std::string, M>& res) const {
1930  // Get default inputs
1931  std::vector<M> res_v(n_out_);
1932  for (casadi_int i=0; i<res_v.size(); ++i) {
1933  res_v[i] = std::numeric_limits<double>::quiet_NaN();
1934  }
1935 
1936  // Assign provided inputs
1937  for (auto&& e : res) {
1938  M a = e.second;
1939  res_v.at(index_out(e.first)) = a;
1940  }
1941  return res_v;
1942  }
1943 
1944  template<typename MatType>
1945  MatType FunctionInternal::ensure_stacked(const MatType& v, const Sparsity& sp, casadi_int n) {
1946  // Check dimensions
1947  if (v.size1() == sp.size1() && v.size2() == n * sp.size2()) {
1948  // Ensure that sparsity is a horizontal multiple of original input, or has no entries
1949  if (v.nnz() != 0 && !v.sparsity().is_stacked(sp, n)) {
1950  return project(v, repmat(sp, 1, n));
1951  }
1952  } else {
1953  // Correct empty sparsity
1954  casadi_assert_dev(v.is_empty());
1955  return MatType(sp.size1(), sp.size2() * n);
1956  }
1957  // No correction needed
1958  return v;
1959  }
1960 
1961 } // namespace casadi
1962 
1964 
1965 #endif // CASADI_FUNCTION_INTERNAL_HPP
casadi::Function::AuxOut
std::map< std::string, std::vector< std::string > > AuxOut
Definition: function.hpp:447
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:1289
casadi
The casadi namespace.
Definition: archiver.hpp:32
casadi::MXVector
std::vector< MX > MXVector
Definition: mx.hpp:1107
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