sparsity_interface.hpp
1 /*
2  * This file is part of CasADi.
3  *
4  * CasADi -- A symbolic framework for dynamic optimization.
5  * Copyright (C) 2010-2023 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_SPARSITY_INTERFACE_HPP
27 #define CASADI_SPARSITY_INTERFACE_HPP
28 
29 #include "casadi_misc.hpp"
30 
31 namespace casadi {
37  struct CASADI_EXPORT SparsityInterfaceCommon {};
38 
39 #ifndef SWIG
50  template<typename MatType>
51  class SparsityInterface : public SparsityInterfaceCommon {
52  protected:
53  // Helper functions
54  inline const MatType& self() const { return static_cast<const MatType&>(*this); }
55  inline MatType& self() { return static_cast<MatType&>(*this); }
56  public:
57 
59  static std::vector< std::vector< MatType > >
60  blocksplit(const MatType& x, const std::vector<casadi_int>& vert_offset,
61  const std::vector<casadi_int>& horz_offset);
62  static std::vector< std::vector< MatType > >
63  blocksplit(const MatType& x, casadi_int vert_incr, casadi_int horz_incr);
64  static MatType veccat(const std::vector< MatType >& x);
65  static MatType vec(const MatType& x);
66  static MatType repmat(const MatType& x, casadi_int n, casadi_int m=1);
67  static std::vector<casadi_int> offset(const std::vector< MatType > &v, bool vert=true);
68  static std::vector< MatType > diagsplit(const MatType& x,
69  const std::vector<casadi_int>& output_offset);
70  static std::vector< MatType > diagsplit(const MatType& x, casadi_int incr);
71  static std::vector< MatType > diagsplit(const MatType& x, casadi_int incr1, casadi_int incr2);
72  static MatType mtimes(const std::vector<MatType> &args);
73  static std::vector<MatType > horzsplit(const MatType& x, casadi_int incr);
74  static std::vector<MatType > vertsplit(const MatType& x, casadi_int incr);
75  static std::vector<MatType > horzsplit_n(const MatType& x, casadi_int n);
76  static std::vector<MatType > vertsplit_n(const MatType& x, casadi_int n);
77  static MatType repmat(const MatType &A, const std::pair<casadi_int, casadi_int>& rc) {
78  return MatType::repmat(A, rc.first, rc.second);
79  }
81 
97  inline friend MatType horzcat(const std::vector<MatType> &v) {
98  return MatType::horzcat(v);
99  }
100 
108  inline friend MatType vertcat(const std::vector<MatType> &v) {
109  return MatType::vertcat(v);
110  }
111 
120  inline friend std::vector<MatType >
121  horzsplit(const MatType &x, const std::vector<casadi_int>& offset) {
122  return MatType::horzsplit(x, offset);
123  }
124 
134  inline friend std::vector<MatType > horzsplit(const MatType& x, casadi_int incr=1) {
135  return MatType::horzsplit(x, incr);
136  }
137 
149  inline friend std::vector<MatType > horzsplit_n(const MatType& x, casadi_int n) {
150  return MatType::horzsplit_n(x, n);
151  }
152 
161  inline friend std::vector<MatType >
162  vertsplit(const MatType& x, const std::vector<casadi_int>& offset) {
163  return MatType::vertsplit(x, offset);
164  }
165 
169  inline friend std::vector<casadi_int > offset(const std::vector<MatType> &v, bool vert=true) {
170  return MatType::offset(v, vert);
171  }
172 
204  inline friend std::vector<MatType > vertsplit(const MatType &x, casadi_int incr=1) {
205  return MatType::vertsplit(x, incr);
206  }
207 
219  inline friend std::vector<MatType > vertsplit_n(const MatType& x, casadi_int n) {
220  return MatType::vertsplit_n(x, n);
221  }
222 
226  inline friend MatType blockcat(const std::vector< std::vector<MatType > > &v) {
227  return MatType::blockcat(v);
228  }
229 
233  inline friend MatType
234  blockcat(const MatType &A, const MatType &B, const MatType &C, const MatType &D) {
235  return vertcat(horzcat(A, B), horzcat(C, D));
236  }
237 
246  inline friend std::vector< std::vector< MatType > >
247  blocksplit(const MatType& x,
248  const std::vector<casadi_int>& vert_offset,
249  const std::vector<casadi_int>& horz_offset) {
250  return MatType::blocksplit(x, vert_offset, horz_offset);
251  }
252 
261  inline friend std::vector< std::vector< MatType > >
262  blocksplit(const MatType& x, casadi_int vert_incr=1, casadi_int horz_incr=1) {
263  return MatType::blocksplit(x, vert_incr, horz_incr);
264  }
265 
269  inline friend MatType diagcat(const std::vector<MatType> &A) {
270  return MatType::diagcat(A);
271  }
272 
283  friend std::vector< MatType >
284  diagsplit(const MatType& x,
285  const std::vector<casadi_int>& output_offset1,
286  const std::vector<casadi_int>& output_offset2) {
287  return MatType::diagsplit(x, output_offset1, output_offset2);
288  }
289 
298  inline friend std::vector< MatType >
299  diagsplit(const MatType& x, const std::vector<casadi_int>& output_offset) {
300  return MatType::diagsplit(x, output_offset);
301  }
302 
310  inline friend std::vector< MatType >
311  diagsplit(const MatType& x, casadi_int incr=1) {
312  return MatType::diagsplit(x, incr);
313  }
314 
323  inline friend std::vector< MatType >
324  diagsplit(const MatType& x, casadi_int incr1, casadi_int incr2) {
325  return MatType::diagsplit(x, incr1, incr2);
326  }
327 
331  inline friend MatType veccat(const std::vector< MatType >& x) {
332  return MatType::veccat(x);
333  }
334 
338  inline friend MatType mtimes(const MatType &x, const MatType &y) {
339  return MatType::mtimes(x, y);
340  }
341 
345  inline friend MatType mtimes(const std::vector<MatType> &args) {
346  return MatType::mtimes(args);
347  }
348 
357  inline friend MatType
358  mac(const MatType &x, const MatType &y, const MatType &z) {
359  return MatType::mac(x, y, z);
360  }
361 
365  inline friend MatType transpose(const MatType& X) {
366  return X.T();
367  }
368 
386  inline friend MatType vec(const MatType& x) {
387  return MatType::vec(x);
388  }
389 
393  inline friend MatType reshape(const MatType& x, casadi_int nrow, casadi_int ncol) {
394  return MatType::reshape(x, nrow, ncol);
395  }
396 
400  inline friend MatType reshape(const MatType& x, std::pair<casadi_int, casadi_int> rc) {
401  return MatType::reshape(x, rc.first, rc.second);
402  }
403 
407  inline friend MatType reshape(const MatType& x, const Sparsity& sp) {
408  return MatType::reshape(x, sp);
409  }
410 
414  inline friend MatType sparsity_cast(const MatType& x, const Sparsity& sp) {
415  return MatType::sparsity_cast(x, sp);
416  }
417 
421  inline friend casadi_int sprank(const MatType& x) {
422  return MatType::sprank(x);
423  }
424 
428  inline friend casadi_int norm_0_mul(const MatType &x, const MatType &y) {
429  return MatType::norm_0_mul(x, y);
430  }
431 
435  inline friend MatType triu(const MatType& x, bool includeDiagonal=true) {
436  return MatType::triu(x, includeDiagonal);
437  }
438 
442  inline friend MatType tril(const MatType& x, bool includeDiagonal=true) {
443  return MatType::tril(x, includeDiagonal);
444  }
445 
451  inline friend MatType kron(const MatType& a, const MatType& b) {
452  return MatType::kron(a, b);
453  }
454 
458  inline friend MatType repmat(const MatType &A, casadi_int n, casadi_int m=1) {
459  return MatType::repmat(A, n, m);
460  }
461 
465  inline friend MatType repmat(const MatType &A, const std::pair<casadi_int, casadi_int>& rc) {
466  return MatType::repmat(A, rc);
467  }
468 
472  inline friend MatType horzcat(const MatType &x, const MatType &y) {
473  return horzcat(std::vector<MatType>{x, y});
474  }
475 
479  inline friend MatType horzcat(const MatType &x, const MatType &y, const MatType &z) {
480  return horzcat(std::vector<MatType>{x, y, z});
481  }
482 
486  inline friend MatType horzcat(const MatType &x, const MatType &y, const MatType &z,
487  const MatType &w) {
488  return horzcat(std::vector<MatType>{x, y, z, w});
489  }
490 
494  inline friend MatType horzcat(const MatType &x, const MatType &y, const MatType &z,
495  const MatType &w, const MatType &v) {
496  return horzcat(std::vector<MatType>{x, y, z, w, v});
497  }
498 
502  inline friend MatType horzcat(const MatType &x, const MatType &y, const MatType &z,
503  const MatType &w, const MatType &v, const MatType &u) {
504  return horzcat(std::vector<MatType>{x, y, z, w, v, u});
505  }
506 
510  inline friend MatType vertcat(const MatType &x, const MatType &y) {
511  return vertcat(std::vector<MatType>{x, y});
512  }
513 
517  inline friend MatType vertcat(const MatType &x, const MatType &y, const MatType &z) {
518  return vertcat(std::vector<MatType>{x, y, z});
519  }
520 
524  inline friend MatType vertcat(const MatType &x, const MatType &y, const MatType &z,
525  const MatType &w) {
526  return vertcat(std::vector<MatType>{x, y, z, w});
527  }
528 
532  inline friend MatType vertcat(const MatType &x, const MatType &y, const MatType &z,
533  const MatType &w, const MatType &v) {
534  return vertcat(std::vector<MatType>{x, y, z, w, v});
535  }
536 
540  inline friend MatType vertcat(const MatType &x, const MatType &y, const MatType &z,
541  const MatType &w, const MatType &v, const MatType &u) {
542  return vertcat(std::vector<MatType>{x, y, z, w, v, u});
543  }
544 
548  inline friend MatType diagcat(const MatType &x, const MatType &y) {
549  return diagcat(std::vector<MatType>{x, y});
550  }
551 
555  inline friend MatType diagcat(const MatType &x, const MatType &y, const MatType &z) {
556  return diagcat(std::vector<MatType>{x, y, z});
557  }
558 
562  inline friend MatType diagcat(const MatType &x, const MatType &y, const MatType &z,
563  const MatType &w) {
564  return diagcat(std::vector<MatType>{x, y, z, w});
565  }
566 
570  inline friend MatType diagcat(const MatType &x, const MatType &y, const MatType &z,
571  const MatType &w, const MatType &v) {
572  return diagcat(std::vector<MatType>{x, y, z, w, v});
573  }
574 
578  inline friend MatType diagcat(const MatType &x, const MatType &y, const MatType &z,
579  const MatType &w, const MatType &v, const MatType &u) {
580  return diagcat(std::vector<MatType>{x, y, z, w, v, u});
581  }
582 
586  inline friend MatType sum1(const MatType &x) { return MatType::sum1(x);}
587 
591  inline friend MatType sum2(const MatType &x) { return MatType::sum2(x);}
592 
594  };
595 #endif // SWIG
596 
597 #ifndef SWIG
598  template<typename MatType>
599  MatType SparsityInterface<MatType>::vec(const MatType& x) {
600  if (x.size2()==1) return x;
601  return reshape(x, x.numel(), 1);
602  }
603 
604  template<typename MatType>
605  MatType SparsityInterface<MatType>::repmat(const MatType& x, casadi_int n, casadi_int m) {
606  if (n==1 && m==1) return x;
607  MatType allrows = vertcat(std::vector<MatType>(n, x));
608  if (n==0) allrows = MatType(0, x.size2());
609  MatType ret = horzcat(std::vector<MatType>(m, allrows));
610  if (m==0) ret = MatType(allrows.size1(), 0);
611  return ret;
612  }
613 
614  template<typename MatType>
615  std::vector< std::vector< MatType > >
616  SparsityInterface<MatType>::blocksplit(const MatType& x,
617  const std::vector<casadi_int>& vert_offset,
618  const std::vector<casadi_int>& horz_offset) {
619  std::vector<MatType> rows = MatType::vertsplit(x, vert_offset);
620  std::vector< std::vector< MatType > > ret;
621  for (auto&& r : rows) ret.push_back(MatType::horzsplit(r, horz_offset));
622  return ret;
623  }
624 
625  template<typename MatType>
626  std::vector< std::vector< MatType > >
627  SparsityInterface<MatType>::blocksplit(const MatType& x,
628  casadi_int vert_incr, casadi_int horz_incr) {
629  casadi_assert_dev(horz_incr>=1);
630  casadi_assert_dev(vert_incr>=1);
631  casadi_int sz1 = x.size1();
632  std::vector<casadi_int> offset1 = range(0, sz1, vert_incr);
633  offset1.push_back(sz1);
634  casadi_int sz2 = x.size2();
635  std::vector<casadi_int> offset2 = range(0, sz2, horz_incr);
636  offset2.push_back(sz2);
637  return blocksplit(x, offset1, offset2);
638  }
639 
640  template<typename MatType>
641  std::vector<casadi_int>
642  SparsityInterface<MatType>::offset(const std::vector< MatType > &v, bool vert) {
643  std::vector<casadi_int> ret(v.size()+1);
644  ret[0]=0;
645  for (casadi_int i=0; i<v.size(); ++i) {
646  ret[i+1] = ret[i] + (vert ? v[i].size1() : v[i].size2());
647  }
648  return ret;
649  }
650 
651  template<typename MatType>
652  MatType SparsityInterface<MatType>::veccat(const std::vector< MatType >& x) {
653  std::vector< MatType > x_vec = x;
654  for (typename std::vector< MatType >::iterator it=x_vec.begin();
655  it!=x_vec.end(); ++it) {
656  *it = vec(*it);
657  }
658  if (x_vec.empty()) {
659  return MatType(0, 1);
660  } else {
661  return vertcat(x_vec);
662  }
663  }
664 
665  template<typename MatType>
666  std::vector< MatType >
667  SparsityInterface<MatType>::diagsplit(const MatType& x,
668  const std::vector<casadi_int>& output_offset) {
669  casadi_assert(x.is_square(), "diagsplit(x,incr)::input must be square but got "
670  + x.dim() + ".");
671  return MatType::diagsplit(x, output_offset, output_offset);
672  }
673 
674  template<typename MatType>
675  std::vector< MatType >
676  SparsityInterface<MatType>::diagsplit(const MatType& x, casadi_int incr) {
677  casadi_assert_dev(incr>=1);
678  casadi_assert(x.is_square(), "diagsplit(x,incr)::input must be square but got "
679  + x.dim() + ".");
680  std::vector<casadi_int> offset2 = range(0, x.size2(), incr);
681  offset2.push_back(x.size2());
682  return MatType::diagsplit(x, offset2);
683  }
684 
685  template<typename MatType>
686  std::vector< MatType >
687  SparsityInterface<MatType>::diagsplit(const MatType& x, casadi_int incr1, casadi_int incr2) {
688  casadi_assert_dev(incr1>=1);
689  casadi_assert_dev(incr2>=1);
690  std::vector<casadi_int> offset1 = range(0, x.size1(), incr1);
691  offset1.push_back(x.size1());
692  std::vector<casadi_int> offset2 = range(0, x.size2(), incr2);
693  offset2.push_back(x.size2());
694  return MatType::diagsplit(x, offset1, offset2);
695  }
696 
697  template<typename MatType>
698  MatType SparsityInterface<MatType>::mtimes(const std::vector<MatType> &args) {
699  casadi_assert(!args.empty(),
700  "mul(std::vector<MatType> &args): "
701  "supplied list must not be empty.");
702  MatType ret = args[0];
703  for (casadi_int i=1; i<args.size(); ++i) ret = MatType::mtimes(ret, args[i]);
704  return ret;
705  }
706 
707  template<typename MatType>
708  std::vector<MatType > SparsityInterface<MatType>::horzsplit(const MatType& x, casadi_int incr) {
709  casadi_assert_dev(incr>=1);
710  casadi_int sz2 = x.size2();
711  std::vector<casadi_int> offset2 = range(0, sz2, incr);
712  offset2.push_back(sz2);
713  return MatType::horzsplit(x, offset2);
714  }
715 
716  template<typename MatType>
717  std::vector<MatType > SparsityInterface<MatType>::vertsplit(const MatType& x, casadi_int incr) {
718  casadi_assert_dev(incr>=1);
719  casadi_int sz1 = x.size1();
720  std::vector<casadi_int> offset1 = range(0, sz1, incr);
721  offset1.push_back(sz1);
722  return MatType::vertsplit(x, offset1);
723  }
724  template<typename MatType>
725  std::vector<MatType > SparsityInterface<MatType>::horzsplit_n(const MatType& x, casadi_int n) {
726  casadi_assert(n>=0, "horzsplit_n(x,n): n (" + str(n) + ") must be non-negative");
727  if (x.size2()==0) return std::vector<MatType>(n, x);
728  casadi_assert(x.size2() % n==0, "horzsplit_n(x,n): x.size2() (" + str(x.size2()) +
729  ") must be a multiple of n (" + str(n) + ")");
730  return horzsplit(x, x.size2()/n);
731  }
732  template<typename MatType>
733  std::vector<MatType > SparsityInterface<MatType>::vertsplit_n(const MatType& x, casadi_int n) {
734  casadi_assert(n>=0, "vertsplit_n(x,n): n (" + str(n) + ") must be non-negative");
735  if (x.size1()==0) return std::vector<MatType>(n, x);
736  casadi_assert(x.size1() % n==0, "vertsplit(x,n): x.size1() (" + str(x.size1()) +
737  ") must be a multiple of n (" + str(n) + ")");
738  return vertsplit(x, x.size1()/n);
739  }
740 
741 #endif // SWIG
742 
743 } // namespace casadi
744 
745 #endif // CASADI_SPARSITY_INTERFACE_HPP
casadi
The casadi namespace.
Definition: binary_mx_impl.hpp:36