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*, std::pair<Function, size_t> >& 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*, std::pair<Function, size_t> >& all_fun,
377  casadi_int max_depth) const;
378 
387  virtual std::vector<bool> which_depends(const std::string& s_in,
388  const std::vector<std::string>& s_out,
389  casadi_int order, bool tr=false) const;
390 
397  virtual Function simplify_passes(
398  const std::vector<std::pair<std::string, casadi_int> >& tasks) const;
399 
401 
404  virtual bool has_spfwd() const { return false;}
405  virtual bool has_sprev() const { return false;}
407 
409 
412  int eval_gen(const double** arg, double** res, casadi_int* iw, double* w, void* mem,
413  bool always_inline, bool never_inline) const;
414  virtual int eval(const double** arg, double** res, casadi_int* iw, double* w, void* mem) const;
416 
420  virtual int eval_sx(const SXElem** arg, SXElem** res,
421  casadi_int* iw, SXElem* w, void* mem, bool always_inline, bool never_inline) const;
422 
426  virtual void eval_mx(const MXVector& arg, MXVector& res,
427  bool always_inline, bool never_inline) const;
428 
430 
433  virtual std::vector<DM> eval_dm(const std::vector<DM>& arg) const;
434  virtual bool has_eval_dm() const { return false;}
436 
438 
441  int eval_gen(const SXElem** arg, SXElem** res, casadi_int* iw, SXElem* w, void* mem,
442  bool always_inline, bool never_inline) const {
443  return eval_sx(arg, res, iw, w, mem, always_inline, never_inline);
444  }
445  int eval_gen(const bvec_t** arg, bvec_t** res, casadi_int* iw, bvec_t* w, void* mem,
446  bool always_inline, bool never_inline) const {
447  return sp_forward(arg, res, iw, w, mem);
448  }
450 
452 
455  void call_gen(const MXVector& arg, MXVector& res, casadi_int npar,
456  bool always_inline, bool never_inline) const;
457 
458  template<typename D>
459  void call_gen(const std::vector<Matrix<D> >& arg, std::vector<Matrix<D> >& res,
460  casadi_int npar, bool always_inline, bool never_inline) const;
462 
466  template<typename M>
467  void call(const std::vector<M>& arg, std::vector<M>& res,
468  bool always_inline, bool never_inline) const;
469 
471 
476  static bool check_mat(const Sparsity& arg, const Sparsity& inp, casadi_int& npar);
478 
480 
488  template<typename M>
489  void check_arg(const std::vector<M>& arg, casadi_int& npar) const;
491 
493 
501  template<typename M>
502  void check_res(const std::vector<M>& res, casadi_int& npar) const;
504 
513  template<typename M> bool
514  matching_arg(const std::vector<M>& arg, casadi_int& npar) const;
515 
524  template<typename M> bool
525  matching_res(const std::vector<M>& arg, casadi_int& npar) const;
526 
530  template<typename M> std::vector<M>
531  replace_arg(const std::vector<M>& arg, casadi_int npar) const;
532 
536  template<typename M> std::vector<M>
537  project_arg(const std::vector<M>& arg, casadi_int npar) const;
538 
542  template<typename M> std::vector<M>
543  project_res(const std::vector<M>& arg, casadi_int npar) const;
544 
548  template<typename M> std::vector<M>
549  replace_res(const std::vector<M>& res, casadi_int npar) const;
550 
554  template<typename M> std::vector<std::vector<M>>
555  replace_fseed(const std::vector<std::vector<M>>& fseed, casadi_int npar) const;
556 
560  template<typename M> std::vector<std::vector<M>>
561  replace_aseed(const std::vector<std::vector<M>>& aseed, casadi_int npar) const;
562 
567  template<typename M>
568  std::map<std::string, M> convert_arg(const std::vector<M>& arg) const;
569  template<typename M>
570  std::vector<M> convert_arg(const std::map<std::string, M>& arg) const;
571  template<typename M>
572  std::map<std::string, M> convert_res(const std::vector<M>& res) const;
573  template<typename M>
574  std::vector<M> convert_res(const std::map<std::string, M>& res) const;
576 
581  std::vector<double> nz_in(const std::vector<DM>& arg) const;
582  std::vector<double> nz_out(const std::vector<DM>& res) const;
583  std::vector<DM> nz_in(const std::vector<double>& arg) const;
584  std::vector<DM> nz_out(const std::vector<double>& res) const;
586 
588 
591  virtual void call_forward(const std::vector<MX>& arg, const std::vector<MX>& res,
592  const std::vector<std::vector<MX> >& fseed,
593  std::vector<std::vector<MX> >& fsens,
594  bool always_inline, bool never_inline) const;
595  virtual void call_forward(const std::vector<SX>& arg, const std::vector<SX>& res,
596  const std::vector<std::vector<SX> >& fseed,
597  std::vector<std::vector<SX> >& fsens,
598  bool always_inline, bool never_inline) const;
600 
602 
605  virtual void call_reverse(const std::vector<MX>& arg, const std::vector<MX>& res,
606  const std::vector<std::vector<MX> >& aseed,
607  std::vector<std::vector<MX> >& asens,
608  bool always_inline, bool never_inline) const;
609  virtual void call_reverse(const std::vector<SX>& arg, const std::vector<SX>& res,
610  const std::vector<std::vector<SX> >& aseed,
611  std::vector<std::vector<SX> >& asens,
612  bool always_inline, bool never_inline) const;
614 
618  std::vector<MX> mapsum_mx(const std::vector<MX > &arg, const std::string& parallelization);
619 
623  virtual bool uses_output() const {return false;}
624 
626 
629  Function jacobian() const;
630  virtual bool has_jacobian() const { return false;}
631  virtual Function get_jacobian(const std::string& name,
632  const std::vector<std::string>& inames,
633  const std::vector<std::string>& onames,
634  const Dict& opts) const;
636 
638 
642  Sparsity& jac_sparsity(casadi_int oind, casadi_int iind, bool compact, bool symmetric) const;
643  virtual bool has_jac_sparsity(casadi_int oind, casadi_int iind) const { return false;}
644  virtual Sparsity get_jac_sparsity(casadi_int oind, casadi_int iind, bool symmetric) const;
646 
648  static std::string forward_name(const std::string& fcn, casadi_int nfwd) {
649  return "fwd" + str(nfwd) + "_" + fcn;
650  }
651 
653  std::string diff_prefix(const std::string& prefix) const;
654 
656 
663  Function forward(casadi_int nfwd) const;
664  virtual bool has_forward(casadi_int nfwd) const { return false;}
665  virtual Function get_forward(casadi_int nfwd, const std::string& name,
666  const std::vector<std::string>& inames,
667  const std::vector<std::string>& onames,
668  const Dict& opts) const;
670 
672  static std::string reverse_name(const std::string& fcn, casadi_int nadj) {
673  return "adj" + str(nadj) + "_" + fcn;
674  }
675 
677 
684  Function reverse(casadi_int nadj) const;
685  virtual bool has_reverse(casadi_int nadj) const { return false;}
686  virtual Function get_reverse(casadi_int nadj, const std::string& name,
687  const std::vector<std::string>& inames,
688  const std::vector<std::string>& onames,
689  const Dict& opts) const;
691 
695  template<typename MatType>
696  static MatType ensure_stacked(const MatType& v, const Sparsity& sp, casadi_int n);
697 
701  virtual Function slice(const std::string& name, const std::vector<casadi_int>& order_in,
702  const std::vector<casadi_int>& order_out, const Dict& opts) const;
703 
707  virtual const Function& oracle() const;
708 
712  bool has_derivative() const;
713 
717  virtual double ad_weight() const;
718 
724  virtual double sp_weight() const;
725 
727 
730  virtual const SX sx_in(casadi_int ind) const;
731  virtual const SX sx_out(casadi_int ind) const;
732  virtual const std::vector<SX> sx_in() const;
733  virtual const std::vector<SX> sx_out() const;
734  virtual const MX mx_in(casadi_int ind) const;
735  virtual const MX mx_out(casadi_int ind) const;
736  virtual const std::vector<MX> mx_in() const;
737  virtual const std::vector<MX> mx_out() const;
738  const DM dm_in(casadi_int ind) const;
739  const DM dm_out(casadi_int ind) const;
740  const std::vector<DM> dm_in() const;
741  const std::vector<DM> dm_out() const;
743 
745  virtual std::vector<MX> free_mx() const;
746 
748  virtual std::vector<SX> free_sx() const;
749 
753  virtual bool has_free() const { return false;}
754 
758  virtual void generate_lifted(Function& vdef_fcn, Function& vinit_fcn) const;
759 
763  virtual casadi_int n_instructions() const;
764 
768  virtual casadi_int instruction_id(casadi_int k) const;
769 
773  virtual std::vector<casadi_int> instruction_input(casadi_int k) const;
774 
778  virtual double instruction_constant(casadi_int k) const;
779 
783  virtual std::vector<casadi_int> instruction_output(casadi_int k) const;
784 
788  virtual casadi_int n_nodes() const;
789 
793  virtual MX instruction_MX(casadi_int k) const;
794 
798  virtual SX instructions_sx() const;
799 
804  Function wrap(const std::string& name) const;
805  Function wrap() const;
807 
812  Function wrap_as_needed(const std::string& name, const Dict& opts) const;
813  Function wrap_as_needed(const Dict& opts) const;
815 
819  Dict cache() const;
820 
824  bool incache(const std::string& fname, Function& f, const std::string& suffix="") const;
825 
829  void tocache(const Function& f, const std::string& suffix="") const;
830 
831 
835  void tocache_if_missing(Function& f, const std::string& suffix="") const;
836 
840  void codegen(CodeGenerator& g, const std::string& fname) const;
841 
845  void codegen_meta(CodeGenerator& g) const;
846 
850  void codegen_sparsities(CodeGenerator& g) const;
851 
855  virtual std::string codegen_name(const CodeGenerator& g, bool ns=true) const;
856 
860  std::string codegen_mem(CodeGenerator& g, const std::string& index="mem") const;
861 
865  virtual void codegen_incref(CodeGenerator& g) const;
866 
870  virtual void codegen_decref(CodeGenerator& g) const;
871 
875  virtual void codegen_alloc_mem(CodeGenerator& g) const;
876 
880  virtual void codegen_init_mem(CodeGenerator& g) const;
881 
885  virtual void codegen_free_mem(CodeGenerator& g) const {}
886 
890  virtual void codegen_checkout(CodeGenerator& g) const;
891 
895  virtual void codegen_release(CodeGenerator& g) const;
896 
900  std::string signature(const std::string& fname) const;
901 
905  std::string signature_unrolled(const std::string& fname) const;
906 
910  virtual void codegen_declarations(CodeGenerator& g) const;
911 
915  virtual void codegen_body(CodeGenerator& g) const;
916 
920  virtual std::string codegen_mem_type() const { return ""; }
921 
925  virtual bool codegen_needs_mem() const { return false; }
926 
932  virtual bool codegen_mem_is_opaque() const { return false; }
933 
937  virtual std::string generate_dependencies(const std::string& fname, const Dict& opts) const;
938 
942  virtual bool has_codegen() const { return false;}
943 
947  virtual void jit_dependencies(const std::string& fname) {}
948 
952  static std::string get_jit_directory(const Dict& jit_options);
953 
957  virtual void export_code(const std::string& lang,
958  std::ostream &stream, const Dict& options) const;
959 
963  void serialize_type(SerializingStream &s) const override;
964 
968  void serialize_body(SerializingStream &s) const override;
969 
973  void disp(std::ostream& stream, bool more) const override;
974 
978  virtual void disp_more(std::ostream& stream) const {}
979 
983  std::string definition() const;
984 
988  void print_dimensions(std::ostream &stream) const;
989 
993  virtual std::vector<std::string> get_free() const;
994 
998  void get_partition(casadi_int iind, casadi_int oind, Sparsity& D1, Sparsity& D2,
999  bool compact, bool symmetric,
1000  bool allow_forward, bool allow_reverse) const;
1001 
1003 
1006  casadi_int nnz_in() const;
1007  casadi_int nnz_in(casadi_int ind) const { return sparsity_in(ind).nnz(); }
1008  casadi_int nnz_out() const;
1009  casadi_int nnz_out(casadi_int ind) const { return sparsity_out(ind).nnz(); }
1011 
1013 
1016  casadi_int numel_in() const;
1017  casadi_int numel_in(casadi_int ind) const { return sparsity_in(ind).numel(); }
1018  casadi_int numel_out(casadi_int ind) const { return sparsity_out(ind).numel(); }
1019  casadi_int numel_out() const;
1021 
1023 
1026  casadi_int size1_in(casadi_int ind) const { return sparsity_in(ind).size1(); }
1027  casadi_int size2_in(casadi_int ind) const { return sparsity_in(ind).size2(); }
1028  casadi_int size1_out(casadi_int ind) const { return sparsity_out(ind).size1(); }
1029  casadi_int size2_out(casadi_int ind) const { return sparsity_out(ind).size2(); }
1030  std::pair<casadi_int, casadi_int> size_in(casadi_int ind) const {
1031  return sparsity_in(ind).size();
1032  }
1033  std::pair<casadi_int, casadi_int> size_out(casadi_int ind) const {
1034  return sparsity_out(ind).size();
1035  }
1037 
1039 
1042  const Sparsity& sparsity_in(casadi_int ind) const { return sparsity_in_.at(ind); }
1043  const Sparsity& sparsity_out(casadi_int ind) const { return sparsity_out_.at(ind); }
1045 
1047 
1050  bool all_scalar() const;
1051 
1053  virtual bool jac_is_symm(casadi_int oind, casadi_int iind) const;
1054 
1056  Sparsity to_compact(casadi_int oind, casadi_int iind, const Sparsity& sp) const;
1057 
1059  Sparsity from_compact(casadi_int oind, casadi_int iind, const Sparsity& sp) const;
1060 
1062  template<bool fwd>
1063  Sparsity get_jac_sparsity_gen(casadi_int oind, casadi_int iind) const;
1064 
1066  Sparsity get_jac_sparsity_hierarchical(casadi_int oind, casadi_int iind) const;
1067 
1071  Sparsity get_jac_sparsity_hierarchical_symm(casadi_int oind, casadi_int iind) const;
1072 
1074  virtual std::vector<MX> symbolic_output(const std::vector<MX>& arg) const;
1075 
1077 
1080  virtual size_t get_n_in();
1081  virtual size_t get_n_out();
1083 
1084 
1086 
1089  virtual std::string get_name_in(casadi_int i);
1090  virtual std::string get_name_out(casadi_int i);
1092 
1096  virtual double get_default_in(casadi_int ind) const {
1097  return 0;
1098  }
1099 
1103  virtual double get_max_in(casadi_int ind) const {
1104  return inf;
1105  }
1106 
1110  virtual double get_min_in(casadi_int ind) const {
1111  return -inf;
1112  }
1113 
1114  virtual std::vector<double> get_nominal_in(casadi_int ind) const {
1115  return std::vector<double>(nnz_in(ind), 1.);
1116  }
1117 
1118  virtual std::vector<double> get_nominal_out(casadi_int ind) const {
1119  return std::vector<double>(nnz_out(ind), 1.);
1120  }
1121 
1125  virtual double get_reltol() const {
1126  return eps;
1127  }
1128 
1132  virtual double get_abstol() const {
1133  return eps;
1134  }
1135 
1139  virtual Sparsity get_sparsity_in(casadi_int i);
1140 
1144  virtual Sparsity get_sparsity_out(casadi_int i);
1145 
1149  virtual bool get_diff_in(casadi_int i) { return true; }
1150 
1154  virtual bool get_diff_out(casadi_int i) { return true; }
1155 
1159  casadi_int index_in(const std::string &name) const {
1160  for (casadi_int i=0; i<name_in_.size(); ++i) {
1161  if (name_in_[i]==name) return i;
1162  }
1163  casadi_error("FunctionInternal::index_in: could not find entry \""
1164  + name + "\". Available names are: " + str(name_in_) + ".");
1165  return -1;
1166  }
1167 
1171  casadi_int index_out(const std::string &name) const {
1172  for (casadi_int i=0; i<name_out_.size(); ++i) {
1173  if (name_out_[i]==name) return i;
1174  }
1175  casadi_error("FunctionInternal::index_out: could not find entry \""
1176  + name + "\". Available names are: " + str(name_out_) + ".");
1177  return -1;
1178  }
1179 
1183  virtual int sp_forward(const bvec_t** arg, bvec_t** res,
1184  casadi_int* iw, bvec_t* w, void* mem) const;
1185 
1193  virtual int eval_activity(const bvec_t** arg, bvec_t** res,
1194  casadi_int* iw, bvec_t* w, void* mem) const;
1195 
1199  virtual int sp_forward_block(const bvec_t** arg, bvec_t** res,
1200  casadi_int* iw, bvec_t* w, void* mem, casadi_int oind, casadi_int iind) const;
1201 
1205  virtual int sp_reverse(bvec_t** arg, bvec_t** res, casadi_int* iw, bvec_t* w, void* mem) const;
1206 
1210  void sz_work(size_t& sz_arg, size_t& sz_res, size_t& sz_iw, size_t& sz_w) const;
1211 
1215  size_t sz_arg() const { return sz_arg_per_ + sz_arg_tmp_;}
1216 
1220  size_t sz_res() const { return sz_res_per_ + sz_res_tmp_;}
1221 
1225  size_t sz_iw() const { return sz_iw_per_ + sz_iw_tmp_;}
1226 
1230  size_t sz_w() const { return sz_w_per_ + sz_w_tmp_;}
1231 
1236  virtual size_t codegen_sz_arg(const CodeGenerator& g) const;
1237  virtual size_t codegen_sz_res(const CodeGenerator& g) const;
1238  virtual size_t codegen_sz_iw(const CodeGenerator& g) const;
1239  virtual size_t codegen_sz_w(const CodeGenerator& g) const;
1241 
1245  void alloc_arg(size_t sz_arg, bool persistent=false);
1246 
1250  void alloc_res(size_t sz_res, bool persistent=false);
1251 
1255  void alloc_iw(size_t sz_iw, bool persistent=false);
1256 
1260  void alloc_w(size_t sz_w, bool persistent=false);
1261 
1265  void alloc(const Function& f, bool persistent=false, int num_threads=1);
1266 
1270  virtual void set_work(void* mem, const double**& arg, double**& res,
1271  casadi_int*& iw, double*& w) const {}
1272 
1276  virtual void set_temp(void* mem, const double** arg, double** res,
1277  casadi_int* iw, double* w) const {}
1278 
1282  void setup(void* mem, const double** arg, double** res, casadi_int* iw, double* w) const;
1283 
1285 
1288  virtual bool fwdViaJac(casadi_int nfwd) const;
1289  virtual bool adjViaJac(casadi_int nadj) const;
1291 
1293  virtual Dict info() const;
1294 
1298  Function map(casadi_int n, const std::string& parallelization) const;
1299 
1303  void generate_in(const std::string& fname, const double** arg) const;
1304  void generate_out(const std::string& fname, double** res) const;
1305 
1306  bool always_inline_, never_inline_;
1307 
1309  size_t n_in_, n_out_;
1310 
1312  std::vector<bool> is_diff_in_, is_diff_out_;
1313 
1315  std::vector<Sparsity> sparsity_in_, sparsity_out_;
1316 
1318  std::vector<std::string> name_in_, name_out_;
1319 
1323  bool jit_;
1324 
1328  bool jit_cleanup_;
1329 
1333  std::string jit_serialize_;
1334 
1338  std::string jit_name_;
1339  std::string jit_directory_;
1340  std::string jit_base_name_;
1341 
1345  bool jit_temp_suffix_;
1346 
1350  eval_t eval_;
1351 
1355  casadi_checkout_t checkout_;
1356 
1360  casadi_release_t release_;
1361 
1365  signal_t incref_, decref_;
1366 
1370  Dict stats_;
1371 
1375  bool has_refcount_;
1376 
1380  bool has_refcount_in_deps_;
1381 
1385  Dict cache_init_;
1386 
1388  mutable WeakCache<std::string, Function> cache_;
1389 
1391  mutable std::vector<Sparsity> jac_sparsity_[2];
1392 
1393 #ifdef CASADI_WITH_THREADSAFE_SYMBOLICS
1395  mutable std::mutex jac_sparsity_mtx_;
1396 #endif // CASADI_WITH_THREADSAFE_SYMBOLICS
1397 
1399  Function derivative_of_;
1400 
1402  void* user_data_;
1403 
1405  std::string compiler_plugin_;
1406  Importer compiler_;
1407  Dict jit_options_;
1408 
1410  double jac_penalty_;
1411 
1412  // Types of derivative calculation permitted
1413  bool enable_forward_, enable_reverse_, enable_jacobian_, enable_fd_;
1414  bool enable_forward_op_, enable_reverse_op_, enable_jacobian_op_, enable_fd_op_;
1415 
1417  double ad_weight_, ad_weight_sp_;
1418 
1420  casadi_int max_num_dir_;
1421 
1423  bool inputs_check_;
1424 
1425  // Finite difference step
1426  Dict fd_options_;
1427 
1428  // Finite difference step size
1429  double fd_step_;
1430 
1431  // Finite difference method
1432  std::string fd_method_;
1433 
1434  // Print input/output
1435  bool print_in_;
1436  bool print_out_;
1437 
1438  // Print canonical form
1439  bool print_canonical_;
1440 
1441  // Warn when number of inputs or outputs exceed this value
1442  casadi_int max_io_;
1443 
1444  // Dump input/output
1445  bool dump_in_, dump_out_, dump_;
1446 
1447  // Directory to dump to
1448  std::string dump_dir_;
1449 
1450  // Format to dump with
1451  std::string dump_format_;
1452 
1453  // Forward/reverse/Jacobian options
1454  Dict forward_options_, reverse_options_, jacobian_options_, der_options_;
1455 
1456  // Store a reference to a custom Jacobian
1457  Function custom_jacobian_;
1458 
1459  // Registered functions
1460  std::vector<Function> registered_functions_;
1461 
1462  // Counter for unique names for dumping inputs and output
1463 #ifdef CASADI_WITH_THREAD
1464  mutable std::atomic<casadi_int> dump_count_;
1465 #else
1466  mutable casadi_int dump_count_;
1467 #endif // CASADI_WITH_THREAD
1468 
1472  virtual bool is_a(const std::string& type, bool recursive) const;
1473 
1477  virtual void merge(const std::vector<MX>& arg,
1478  std::vector<MX>& subs_from, std::vector<MX>& subs_to) const;
1479 
1483  template<typename MatType>
1484  static bool purgable(const std::vector<MatType>& seed);
1485 
1489  template<typename MatType>
1490  std::vector<std::vector<MatType> >
1491  fwd_seed(casadi_int nfwd) const;
1492 
1496  template<typename MatType>
1497  std::vector<std::vector<MatType> >
1498  symbolicAdjSeed(casadi_int nadj, const std::vector<MatType>& v) const;
1499 
1500  static std::string string_from_UnifiedReturnStatus(UnifiedReturnStatus status);
1501 
1505  explicit FunctionInternal(DeserializingStream& e);
1506 
1510  static Function deserialize(DeserializingStream& s);
1511  static std::map<std::string, ProtoFunction* (*)(DeserializingStream&)> deserialize_map;
1512 
1516  void print_in(std::ostream &stream, const double** arg, bool truncate) const;
1517 
1521  void print_out(std::ostream &stream, double** res, bool truncate) const;
1522 
1526  static void print_canonical(std::ostream &stream, const Sparsity& sp, const double* nz);
1527 
1531  static void print_canonical(std::ostream &stream, casadi_int sz, const double* nz);
1532 
1536  static void print_canonical(std::ostream &stream, double a);
1537 
1538  protected:
1542  void set_jac_sparsity(casadi_int oind, casadi_int iind, const Sparsity& sp);
1543 
1544  private:
1545  // @{
1547  casadi_int get_dump_id() const;
1548  void dump_in(casadi_int id, const double** arg) const;
1549  void dump_out(casadi_int id, double** res) const;
1550  void dump() const;
1551  // @}
1552 
1556  size_t sz_arg_per_, sz_res_per_, sz_iw_per_, sz_w_per_;
1557 
1561  size_t sz_arg_tmp_, sz_res_tmp_, sz_iw_tmp_, sz_w_tmp_;
1562  };
1563 
1564  // Template implementations
1565  template<typename MatType>
1566  bool FunctionInternal::purgable(const std::vector<MatType>& v) {
1567  for (auto i=v.begin(); i!=v.end(); ++i) {
1568  if (!i->is_zero()) return false;
1569  }
1570  return true;
1571  }
1572 
1573  template<typename MatType>
1574  std::vector<std::vector<MatType> >
1575  FunctionInternal::
1576  fwd_seed(casadi_int nfwd) const {
1577  std::vector<std::vector<MatType>> fseed(nfwd);
1578  for (casadi_int dir=0; dir<nfwd; ++dir) {
1579  fseed[dir].resize(n_in_);
1580  for (casadi_int iind=0; iind<n_in_; ++iind) {
1581  std::string n = "f" + str(dir) + "_" + name_in_[iind];
1582  Sparsity sp = is_diff_in_[iind] ? sparsity_in(iind) : Sparsity(size_in(iind));
1583  fseed[dir][iind] = MatType::sym(n, sp);
1584  }
1585  }
1586  return fseed;
1587  }
1588 
1589  template<typename MatType>
1590  std::vector<std::vector<MatType> >
1591  FunctionInternal::
1592  symbolicAdjSeed(casadi_int nadj, const std::vector<MatType>& v) const {
1593  std::vector<std::vector<MatType> > aseed(nadj, v);
1594  for (casadi_int dir=0; dir<nadj; ++dir) {
1595  // Replace symbolic inputs
1596  casadi_int oind=0;
1597  for (typename std::vector<MatType>::iterator i=aseed[dir].begin();
1598  i!=aseed[dir].end();
1599  ++i, ++oind) {
1600  // Name of the adjoint seed
1601  std::stringstream ss;
1602  ss << "a";
1603  if (nadj>1) ss << dir << "_";
1604  ss << oind;
1605 
1606  // Save to matrix
1607  *i = MatType::sym(ss.str(), is_diff_out_[oind] ? i->sparsity() : Sparsity(i->size()));
1608 
1609  }
1610  }
1611  return aseed;
1612  }
1613 
1614  template<typename M>
1615  void FunctionInternal::call(const std::vector<M>& arg, std::vector<M>& res,
1616  bool always_inline, bool never_inline) const {
1617  // If all inputs are scalar ...
1618  if (all_scalar()) {
1619  // ... and some arguments are matrix-valued with matching dimensions ...
1620  bool matrix_call = false;
1621  std::pair<casadi_int, casadi_int> sz;
1622  for (auto&& a : arg) {
1623  if (!a.is_scalar() && !a.is_empty()) {
1624  if (!matrix_call) {
1625  // Matrix call
1626  matrix_call = true;
1627  sz = a.size();
1628  } else if (a.size()!=sz) {
1629  // Not same dimensions
1630  matrix_call = false;
1631  break;
1632  }
1633  }
1634  }
1635 
1636  // ... then, call multiple times
1637  if (matrix_call) {
1638  // Start with zeros
1639  res.resize(n_out_);
1640  M z = M::zeros(sz);
1641  for (auto&& a : res) a = z;
1642  // Call multiple times
1643  std::vector<M> arg1 = arg, res1;
1644  for (casadi_int c=0; c<sz.second; ++c) {
1645  for (casadi_int r=0; r<sz.first; ++r) {
1646  // Get scalar arguments
1647  for (casadi_int i=0; i<arg.size(); ++i) {
1648  if (arg[i].size()==sz) arg1[i] = arg[i](r, c);
1649  }
1650  // Call recursively with scalar arguments
1651  call(arg1, res1, always_inline, never_inline);
1652  // Get results
1653  casadi_assert_dev(res.size() == res1.size());
1654  for (casadi_int i=0; i<res.size(); ++i) res[i](r, c) = res1[i];
1655  }
1656  }
1657  // All elements assigned
1658  return;
1659  }
1660  }
1661 
1662  // Check if inputs need to be replaced
1663  casadi_int npar = 1;
1664  if (!matching_arg(arg, npar)) {
1665  return call(replace_arg(arg, npar), res, always_inline, never_inline);
1666  }
1667 
1668  // Call the type-specific method
1669  call_gen(arg, res, npar, always_inline, never_inline);
1670  }
1671 
1672  template<typename M>
1673  std::vector<M> FunctionInternal::
1674  project_arg(const std::vector<M>& arg, casadi_int npar) const {
1675  casadi_assert_dev(arg.size()==n_in_);
1676 
1677  // Which arguments require mapped evaluation
1678  std::vector<bool> mapped(n_in_);
1679  for (casadi_int i=0; i<n_in_; ++i) {
1680  mapped[i] = arg[i].size2()!=size2_in(i);
1681  }
1682 
1683  // Check if matching input sparsity
1684  std::vector<bool> matching(n_in_);
1685  bool any_mismatch = false;
1686  for (casadi_int i=0; i<n_in_; ++i) {
1687  if (mapped[i]) {
1688  matching[i] = arg[i].sparsity().is_stacked(sparsity_in(i), npar);
1689  } else {
1690  matching[i] = arg[i].sparsity()==sparsity_in(i);
1691  }
1692  any_mismatch = any_mismatch || !matching[i];
1693  }
1694 
1695  // Correct input sparsity
1696  if (any_mismatch) {
1697  std::vector<M> arg2(arg);
1698  for (casadi_int i=0; i<n_in_; ++i) {
1699  if (!matching[i]) {
1700  if (mapped[i]) {
1701  arg2[i] = project(arg2[i], repmat(sparsity_in(i), 1, npar));
1702  } else {
1703  arg2[i] = project(arg2[i], sparsity_in(i));
1704  }
1705  }
1706  }
1707  return arg2;
1708  }
1709  return arg;
1710  }
1711 
1712  template<typename M>
1713  std::vector<M> FunctionInternal::
1714  project_res(const std::vector<M>& arg, casadi_int npar) const {
1715  return arg;
1716  }
1717 
1718  template<typename D>
1719  void FunctionInternal::
1720  call_gen(const std::vector<Matrix<D> >& arg, std::vector<Matrix<D> >& res,
1721  casadi_int npar, bool always_inline, bool never_inline) const {
1722  std::vector< Matrix<D> > arg2 = project_arg(arg, npar);
1723 
1724  // Which arguments require mapped evaluation
1725  std::vector<bool> mapped(n_in_);
1726  for (casadi_int i=0; i<n_in_; ++i) {
1727  mapped[i] = arg[i].size2()!=size2_in(i);
1728  }
1729 
1730  // Allocate results
1731  res.resize(n_out_);
1732  for (casadi_int i=0; i<n_out_; ++i) {
1733  if (!res[i].sparsity().is_stacked(sparsity_out(i), npar)) {
1734  res[i] = Matrix<D>::zeros(repmat(sparsity_out(i), 1, npar));
1735  }
1736  }
1737 
1738  // Allocate temporary memory if needed
1739  std::vector<casadi_int> iw_tmp(sz_iw());
1740  std::vector<D> w_tmp(sz_w());
1741 
1742  // Get pointers to input arguments
1743  std::vector<const D*> argp(sz_arg());
1744  for (casadi_int i=0; i<n_in_; ++i) argp[i]=get_ptr(arg2[i]);
1745 
1746  // Get pointers to output arguments
1747  std::vector<D*> resp(sz_res());
1748  for (casadi_int i=0; i<n_out_; ++i) resp[i]=get_ptr(res[i]);
1749 
1750  // For all parallel calls
1751  for (casadi_int p=0; p<npar; ++p) {
1752  // Call memory-less
1753  if (eval_gen(get_ptr(argp), get_ptr(resp),
1754  get_ptr(iw_tmp), get_ptr(w_tmp), memory(0),
1755  always_inline, never_inline)) {
1756  if (error_on_fail_) casadi_error("Evaluation failed");
1757  }
1758  // Update offsets
1759  if (p==npar-1) break;
1760  for (casadi_int i=0; i<n_in_; ++i) if (mapped[i]) argp[i] += nnz_in(i);
1761  for (casadi_int i=0; i<n_out_; ++i) resp[i] += nnz_out(i);
1762  }
1763  }
1764 
1765  template<typename M>
1766  void FunctionInternal::check_arg(const std::vector<M>& arg, casadi_int& npar) const {
1767  casadi_assert(arg.size()==n_in_, "Incorrect number of inputs: Expected "
1768  + str(n_in_) + ", got " + str(arg.size()));
1769  for (casadi_int i=0; i<n_in_; ++i) {
1770  if (!check_mat(arg[i].sparsity(), sparsity_in(i), npar)) {
1771  // Dimensions
1772  std::string d_arg = str(arg[i].size1()) + "-by-" + str(arg[i].size2());
1773  std::string d_in = str(size1_in(i)) + "-by-" + str(size2_in(i));
1774  std::string e = "Input " + str(i) + " (" + name_in_[i] + ") has mismatching shape. "
1775  "Got " + d_arg + ". Allowed dimensions, in general, are:\n"
1776  " - The input dimension N-by-M (here " + d_in + ")\n"
1777  " - A scalar, i.e. 1-by-1\n"
1778  " - M-by-N if N=1 or M=1 (i.e. a transposed vector)\n"
1779  " - N-by-M1 if K*M1=M for some K (argument repeated horizontally)\n";
1780  if (npar!=-1) {
1781  e += " - N-by-P*M, indicating evaluation with multiple arguments (P must be a "
1782  "multiple of " + str(npar) + " for consistency with previous inputs)";
1783  }
1784  casadi_error(e);
1785  }
1786  }
1787  }
1788 
1789  template<typename M>
1790  void FunctionInternal::check_res(const std::vector<M>& res, casadi_int& npar) const {
1791  casadi_assert(res.size()==n_out_, "Incorrect number of outputs: Expected "
1792  + str(n_out_) + ", got " + str(res.size()));
1793  for (casadi_int i=0; i<n_out_; ++i) {
1794  casadi_assert(check_mat(res[i].sparsity(), sparsity_out(i), npar),
1795  "Output " + str(i) + " (" + name_out_[i] + ") has mismatching shape. "
1796  "Expected " + str(size_out(i)) + ", got " + str(res[i].size()));
1797  }
1798  }
1799 
1800  template<typename M>
1801  bool FunctionInternal::matching_arg(const std::vector<M>& arg, casadi_int& npar) const {
1802  check_arg(arg, npar);
1803  for (casadi_int i=0; i<n_in_; ++i) {
1804  if (arg.at(i).size1()!=size1_in(i)) return false;
1805  if (arg.at(i).size2()!=size2_in(i) && arg.at(i).size2()!=npar*size2_in(i)) return false;
1806  }
1807  return true;
1808  }
1809 
1810  template<typename M>
1811  bool FunctionInternal::matching_res(const std::vector<M>& res, casadi_int& npar) const {
1812  check_res(res, npar);
1813  for (casadi_int i=0; i<n_out_; ++i) {
1814  if (res.at(i).size1()!=size1_out(i)) return false;
1815  if (res.at(i).size2()!=size2_out(i) && res.at(i).size2()!=npar*size2_out(i)) return false;
1816  }
1817  return true;
1818  }
1819 
1820  template<typename M>
1821  M replace_mat(const M& arg, const Sparsity& inp, casadi_int npar) {
1822  if (arg.size()==inp.size()) {
1823  // Matching dimensions already
1824  return arg;
1825  } else if (arg.is_scalar()) {
1826  // Scalar assign means set all
1827  return M(inp, arg);
1828  } else if (arg.is_vector() && inp.size()==std::make_pair(arg.size2(), arg.size1())) {
1829  // Transpose vector
1830  return arg.T();
1831  } else if (arg.size1()==inp.size1() && arg.size2()>0 && inp.size2()>0
1832  && inp.size2()%arg.size2()==0) {
1833  // Horizontal repmat
1834  return repmat(arg, 1, inp.size2()/arg.size2());
1835  } else if (npar!=-1 && arg.size1()==inp.size1() && arg.size2()>0 && inp.size2()>0
1836  && (npar*inp.size2())%arg.size2()==0) {
1837  // Multiple evaluation: grow argument horizontally to npar*inp columns
1838  return repmat(arg, 1, (npar*inp.size2())/arg.size2());
1839  } else {
1840  // Empty matrix means set zero (kept last so a 0-by-N argument is first given the
1841  // chance to be recognised as a parallel/repmat call above)
1842  casadi_assert_dev(arg.is_empty());
1843  return M(inp.size());
1844  }
1845  }
1846 
1847  template<typename M>
1848  std::vector<M> FunctionInternal::
1849  replace_arg(const std::vector<M>& arg, casadi_int npar) const {
1850  std::vector<M> r(arg.size());
1851  for (casadi_int i=0; i<r.size(); ++i) r[i] = replace_mat(arg[i], sparsity_in(i), npar);
1852  return r;
1853  }
1854 
1855  template<typename M>
1856  std::vector<M> FunctionInternal::
1857  replace_res(const std::vector<M>& res, casadi_int npar) const {
1858  std::vector<M> r(res.size());
1859  for (casadi_int i=0; i<r.size(); ++i) r[i] = replace_mat(res[i], sparsity_out(i), npar);
1860  return r;
1861  }
1862 
1863  template<typename M>
1864  std::vector<std::vector<M> > FunctionInternal::
1865  replace_fseed(const std::vector<std::vector<M> >& fseed, casadi_int npar) const {
1866  std::vector<std::vector<M> > r(fseed.size());
1867  for (casadi_int d=0; d<r.size(); ++d) r[d] = replace_arg(fseed[d], npar);
1868  return r;
1869  }
1870 
1871  template<typename M>
1872  std::vector<std::vector<M> > FunctionInternal::
1873  replace_aseed(const std::vector<std::vector<M> >& aseed, casadi_int npar) const {
1874  std::vector<std::vector<M> > r(aseed.size());
1875  for (casadi_int d=0; d<r.size(); ++d) r[d] = replace_res(aseed[d], npar);
1876  return r;
1877  }
1878 
1879  template<typename M>
1880  std::map<std::string, M> FunctionInternal::
1881  convert_arg(const std::vector<M>& arg) const {
1882  casadi_assert(arg.size()==n_in_, "Incorrect number of inputs: Expected "
1883  + str(n_in_) + ", got " + str(arg.size()));
1884  std::map<std::string, M> ret;
1885  for (casadi_int i=0;i<n_in_;++i) {
1886  ret[name_in_[i]] = arg[i];
1887  }
1888  return ret;
1889  }
1890 
1891  template<typename M>
1892  std::vector<M> FunctionInternal::
1893  convert_arg(const std::map<std::string, M>& arg) const {
1894  // Get default inputs
1895  std::vector<M> arg_v(n_in_);
1896  for (casadi_int i=0; i<arg_v.size(); ++i) {
1897  arg_v[i] = get_default_in(i);
1898  }
1899 
1900  // Assign provided inputs
1901  for (auto&& e : arg) {
1902  arg_v.at(index_in(e.first)) = e.second;
1903  }
1904 
1905  return arg_v;
1906  }
1907 
1908  template<typename M>
1909  std::map<std::string, M> FunctionInternal::
1910  convert_res(const std::vector<M>& res) const {
1911  casadi_assert(res.size()==n_out_, "Incorrect number of outputs: Expected "
1912  + str(n_out_) + ", got " + str(res.size()));
1913  std::map<std::string, M> ret;
1914  for (casadi_int i=0;i<n_out_;++i) {
1915  ret[name_out_[i]] = res[i];
1916  }
1917  return ret;
1918  }
1919 
1920  template<typename M>
1921  std::vector<M> FunctionInternal::
1922  convert_res(const std::map<std::string, M>& res) const {
1923  // Get default inputs
1924  std::vector<M> res_v(n_out_);
1925  for (casadi_int i=0; i<res_v.size(); ++i) {
1926  res_v[i] = std::numeric_limits<double>::quiet_NaN();
1927  }
1928 
1929  // Assign provided inputs
1930  for (auto&& e : res) {
1931  M a = e.second;
1932  res_v.at(index_out(e.first)) = a;
1933  }
1934  return res_v;
1935  }
1936 
1937  template<typename MatType>
1938  MatType FunctionInternal::ensure_stacked(const MatType& v, const Sparsity& sp, casadi_int n) {
1939  // Check dimensions
1940  if (v.size1() == sp.size1() && v.size2() == n * sp.size2()) {
1941  // Ensure that sparsity is a horizontal multiple of original input, or has no entries
1942  if (v.nnz() != 0 && !v.sparsity().is_stacked(sp, n)) {
1943  return project(v, repmat(sp, 1, n));
1944  }
1945  } else {
1946  // Correct empty sparsity
1947  casadi_assert_dev(v.is_empty());
1948  return MatType(sp.size1(), sp.size2() * n);
1949  }
1950  // No correction needed
1951  return v;
1952  }
1953 
1954 } // namespace casadi
1955 
1957 
1958 #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