Sleipnir C++ API
Loading...
Searching...
No Matches
variable_matrix.hpp
1// Copyright (c) Sleipnir contributors
2
3#pragma once
4
5#include <algorithm>
6#include <concepts>
7#include <initializer_list>
8#include <iterator>
9#include <span>
10#include <utility>
11#include <vector>
12
13#include <Eigen/Core>
14#include <Eigen/QR>
15#include <gch/small_vector.hpp>
16
17#include "sleipnir/autodiff/sleipnir_base.hpp"
18#include "sleipnir/autodiff/slice.hpp"
19#include "sleipnir/autodiff/variable.hpp"
20#include "sleipnir/autodiff/variable_block.hpp"
21#include "sleipnir/util/assert.hpp"
22#include "sleipnir/util/concepts.hpp"
23#include "sleipnir/util/empty.hpp"
24#include "sleipnir/util/function_ref.hpp"
25#include "sleipnir/util/symbol_exports.hpp"
26
27namespace slp {
28
34template <typename Scalar_>
36 public:
40 using Scalar = Scalar_;
41
45 VariableMatrix() = default;
46
53 explicit VariableMatrix(int rows) : VariableMatrix{rows, 1} {}
54
61 VariableMatrix(int rows, int cols) : m_rows{rows}, m_cols{cols} {
62 m_storage.reserve(rows * cols);
63 for (int index = 0; index < rows * cols; ++index) {
64 m_storage.emplace_back();
65 }
66 }
67
75 : m_rows{rows}, m_cols{cols} {
76 m_storage.reserve(rows * cols);
77 for (int index = 0; index < rows * cols; ++index) {
78 m_storage.emplace_back(nullptr);
79 }
80 }
81
88 std::initializer_list<std::initializer_list<Variable<Scalar>>> list) {
89 // Get row and column counts for destination matrix
90 m_rows = list.size();
91 m_cols = 0;
92 if (list.size() > 0) {
93 m_cols = list.begin()->size();
94 }
95
96 // Assert all column counts are the same
97 for ([[maybe_unused]]
98 const auto& row : list) {
99 slp_assert(static_cast<int>(row.size()) == m_cols);
100 }
101
102 m_storage.reserve(rows() * cols());
103 for (const auto& row : list) {
104 std::ranges::copy(row, std::back_inserter(m_storage));
105 }
106 }
107
115 // NOLINTNEXTLINE (google-explicit-constructor)
116 VariableMatrix(const std::vector<std::vector<Scalar>>& list) {
117 // Get row and column counts for destination matrix
118 m_rows = list.size();
119 m_cols = 0;
120 if (list.size() > 0) {
121 m_cols = list.begin()->size();
122 }
123
124 // Assert all column counts are the same
125 for ([[maybe_unused]]
126 const auto& row : list) {
127 slp_assert(static_cast<int>(row.size()) == m_cols);
128 }
129
130 m_storage.reserve(rows() * cols());
131 for (const auto& row : list) {
132 std::ranges::copy(row, std::back_inserter(m_storage));
133 }
134 }
135
143 // NOLINTNEXTLINE (google-explicit-constructor)
144 VariableMatrix(const std::vector<std::vector<Variable<Scalar>>>& list) {
145 // Get row and column counts for destination matrix
146 m_rows = list.size();
147 m_cols = 0;
148 if (list.size() > 0) {
149 m_cols = list.begin()->size();
150 }
151
152 // Assert all column counts are the same
153 for ([[maybe_unused]]
154 const auto& row : list) {
155 slp_assert(static_cast<int>(row.size()) == m_cols);
156 }
157
158 m_storage.reserve(rows() * cols());
159 for (const auto& row : list) {
160 std::ranges::copy(row, std::back_inserter(m_storage));
161 }
162 }
163
169 template <typename Derived>
170 // NOLINTNEXTLINE (google-explicit-constructor)
171 VariableMatrix(const Eigen::MatrixBase<Derived>& values)
172 : m_rows{static_cast<int>(values.rows())},
173 m_cols{static_cast<int>(values.cols())} {
174 m_storage.reserve(values.rows() * values.cols());
175 for (int row = 0; row < values.rows(); ++row) {
176 for (int col = 0; col < values.cols(); ++col) {
177 m_storage.emplace_back(values(row, col));
178 }
179 }
180 }
181
187 template <typename Derived>
188 // NOLINTNEXTLINE (google-explicit-constructor)
189 VariableMatrix(const Eigen::DiagonalBase<Derived>& values)
190 : m_rows{static_cast<int>(values.rows())},
191 m_cols{static_cast<int>(values.cols())} {
192 m_storage.reserve(values.rows() * values.cols());
193 for (int row = 0; row < values.rows(); ++row) {
194 for (int col = 0; col < values.cols(); ++col) {
195 if (row == col) {
196 m_storage.emplace_back(values.diagonal()[row]);
197 } else {
198 m_storage.emplace_back(Scalar(0));
199 }
200 }
201 }
202 }
203
209 // NOLINTNEXTLINE (google-explicit-constructor)
210 VariableMatrix(const Variable<Scalar>& variable) : m_rows{1}, m_cols{1} {
211 m_storage.emplace_back(variable);
212 }
213
219 // NOLINTNEXTLINE (google-explicit-constructor)
220 VariableMatrix(Variable<Scalar>&& variable) : m_rows{1}, m_cols{1} {
221 m_storage.emplace_back(std::move(variable));
222 }
223
229 // NOLINTNEXTLINE (google-explicit-constructor)
231 : m_rows{values.rows()}, m_cols{values.cols()} {
232 m_storage.reserve(rows() * cols());
233 for (int row = 0; row < rows(); ++row) {
234 for (int col = 0; col < cols(); ++col) {
235 m_storage.emplace_back(values[row, col]);
236 }
237 }
238 }
239
245 // NOLINTNEXTLINE (google-explicit-constructor)
247 : m_rows{values.rows()}, m_cols{values.cols()} {
248 m_storage.reserve(rows() * cols());
249 for (int row = 0; row < rows(); ++row) {
250 for (int col = 0; col < cols(); ++col) {
251 m_storage.emplace_back(values[row, col]);
252 }
253 }
254 }
255
261 explicit VariableMatrix(std::span<const Variable<Scalar>> values)
262 : m_rows{static_cast<int>(values.size())}, m_cols{1} {
263 m_storage.reserve(rows() * cols());
264 for (int row = 0; row < rows(); ++row) {
265 for (int col = 0; col < cols(); ++col) {
266 m_storage.emplace_back(values[row * cols() + col]);
267 }
268 }
269 }
270
278 VariableMatrix(std::span<const Variable<Scalar>> values, int rows, int cols)
279 : m_rows{rows}, m_cols{cols} {
280 slp_assert(static_cast<int>(values.size()) == rows * cols);
281 m_storage.reserve(rows * cols);
282 for (int row = 0; row < rows; ++row) {
283 for (int col = 0; col < cols; ++col) {
284 m_storage.emplace_back(values[row * cols + col]);
285 }
286 }
287 }
288
295 template <typename Derived>
296 VariableMatrix& operator=(const Eigen::MatrixBase<Derived>& values) {
297 slp_assert(rows() == values.rows() && cols() == values.cols());
298
299 for (int row = 0; row < values.rows(); ++row) {
300 for (int col = 0; col < values.cols(); ++col) {
301 (*this)[row, col] = values(row, col);
302 }
303 }
304
305 return *this;
306 }
307
317 slp_assert(rows() == 1 && cols() == 1);
318
319 (*this)[0, 0] = value;
320
321 return *this;
322 }
323
329 template <typename Derived>
330 requires std::same_as<typename Derived::Scalar, Scalar>
331 void set_value(const Eigen::MatrixBase<Derived>& values) {
332 slp_assert(rows() == values.rows() && cols() == values.cols());
333
334 for (int row = 0; row < values.rows(); ++row) {
335 for (int col = 0; col < values.cols(); ++col) {
336 (*this)[row, col].set_value(values(row, col));
337 }
338 }
339 }
340
349 slp_assert(row >= 0 && row < rows());
350 slp_assert(col >= 0 && col < cols());
351 return m_storage[row * cols() + col];
352 }
353
361 const Variable<Scalar>& operator[](int row, int col) const {
362 slp_assert(row >= 0 && row < rows());
363 slp_assert(col >= 0 && col < cols());
364 return m_storage[row * cols() + col];
365 }
366
374 slp_assert(index >= 0 && index < rows() * cols());
375 return m_storage[index];
376 }
377
384 const Variable<Scalar>& operator[](int index) const {
385 slp_assert(index >= 0 && index < rows() * cols());
386 return m_storage[index];
387 }
388
399 int block_rows, int block_cols) {
400 slp_assert(row_offset >= 0 && row_offset <= rows());
401 slp_assert(col_offset >= 0 && col_offset <= cols());
402 slp_assert(block_rows >= 0 && block_rows <= rows() - row_offset);
403 slp_assert(block_cols >= 0 && block_cols <= cols() - col_offset);
405 }
406
417 int col_offset,
418 int block_rows,
419 int block_cols) const {
420 slp_assert(row_offset >= 0 && row_offset <= rows());
421 slp_assert(col_offset >= 0 && col_offset <= cols());
422 slp_assert(block_rows >= 0 && block_rows <= rows() - row_offset);
423 slp_assert(block_cols >= 0 && block_cols <= cols() - col_offset);
425 }
426
440
449 Slice col_slice) const {
450 int row_slice_length = row_slice.adjust(rows());
451 int col_slice_length = col_slice.adjust(cols());
452 return VariableBlock{*this, std::move(row_slice), row_slice_length,
453 std::move(col_slice), col_slice_length};
454 }
455
476
495
504 slp_assert(cols() == 1);
505 slp_assert(offset >= 0 && offset < rows());
506 slp_assert(length >= 0 && length <= rows() - offset);
507 return block(offset, 0, length, 1);
508 }
509
518 int length) const {
519 slp_assert(cols() == 1);
520 slp_assert(offset >= 0 && offset < rows());
521 slp_assert(length >= 0 && length <= rows() - offset);
522 return block(offset, 0, length, 1);
523 }
524
532 slp_assert(row >= 0 && row < rows());
533 return block(row, 0, 1, cols());
534 }
535
543 slp_assert(row >= 0 && row < rows());
544 return block(row, 0, 1, cols());
545 }
546
554 slp_assert(col >= 0 && col < cols());
555 return block(0, col, rows(), 1);
556 }
557
565 slp_assert(col >= 0 && col < cols());
566 return block(0, col, rows(), 1);
567 }
568
575 template <EigenMatrixLike LHS, SleipnirMatrixLike<Scalar> RHS>
577 slp_assert(lhs.cols() == rhs.rows());
578
579 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), rhs.cols());
580
581 for (int i = 0; i < lhs.rows(); ++i) {
582 for (int j = 0; j < rhs.cols(); ++j) {
583 Variable sum{Scalar(0)};
584 for (int k = 0; k < lhs.cols(); ++k) {
585 sum += lhs(i, k) * rhs[k, j];
586 }
587 result[i, j] = sum;
588 }
589 }
590
591 return result;
592 }
593
600 template <SleipnirMatrixLike<Scalar> LHS, EigenMatrixLike RHS>
602 slp_assert(lhs.cols() == rhs.rows());
603
604 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), rhs.cols());
605
606 for (int i = 0; i < lhs.rows(); ++i) {
607 for (int j = 0; j < rhs.cols(); ++j) {
608 Variable sum{Scalar(0)};
609 for (int k = 0; k < lhs.cols(); ++k) {
610 sum += lhs[i, k] * rhs(k, j);
611 }
612 result[i, j] = sum;
613 }
614 }
615
616 return result;
617 }
618
625 template <SleipnirMatrixLike<Scalar> LHS, SleipnirMatrixLike<Scalar> RHS>
627 slp_assert(lhs.cols() == rhs.rows());
628
629 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), rhs.cols());
630
631 for (int i = 0; i < lhs.rows(); ++i) {
632 for (int j = 0; j < rhs.cols(); ++j) {
633 Variable sum{Scalar(0)};
634 for (int k = 0; k < lhs.cols(); ++k) {
635 sum += lhs[i, k] * rhs[k, j];
636 }
637 result[i, j] = sum;
638 }
639 }
640
641 return result;
642 }
643
650 template <EigenMatrixLike LHS>
652 const Variable<Scalar>& rhs) {
653 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
654
655 for (int row = 0; row < result.rows(); ++row) {
656 for (int col = 0; col < result.cols(); ++col) {
657 result[row, col] = lhs(row, col) * rhs;
658 }
659 }
660
661 return result;
662 }
663
670 template <SleipnirMatrixLike<Scalar> LHS, ScalarLike RHS>
672 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
673
674 for (int row = 0; row < result.rows(); ++row) {
675 for (int col = 0; col < result.cols(); ++col) {
676 result[row, col] = lhs[row, col] * rhs;
677 }
678 }
679
680 return result;
681 }
682
689 template <EigenMatrixLike RHS>
691 const RHS& rhs) {
692 VariableMatrix<Scalar> result(detail::empty, rhs.rows(), rhs.cols());
693
694 for (int row = 0; row < result.rows(); ++row) {
695 for (int col = 0; col < result.cols(); ++col) {
696 result[row, col] = rhs(row, col) * lhs;
697 }
698 }
699
700 return result;
701 }
702
709 template <ScalarLike LHS, SleipnirMatrixLike<Scalar> RHS>
711 VariableMatrix<Scalar> result(detail::empty, rhs.rows(), rhs.cols());
712
713 for (int row = 0; row < result.rows(); ++row) {
714 for (int col = 0; col < result.cols(); ++col) {
715 result[row, col] = rhs[row, col] * lhs;
716 }
717 }
718
719 return result;
720 }
721
729 slp_assert(cols() == rhs.rows() && cols() == rhs.cols());
730
731 for (int i = 0; i < rows(); ++i) {
732 for (int j = 0; j < rhs.cols(); ++j) {
733 Variable sum{Scalar(0)};
734 for (int k = 0; k < cols(); ++k) {
735 if constexpr (SleipnirMatrixLike<decltype(rhs), Scalar>) {
736 sum += (*this)[i, k] * rhs[k, j];
737 } else {
738 sum += (*this)[i, k] * rhs(k, j);
739 }
740 }
741 (*this)[i, j] = sum;
742 }
743 }
744
745 return *this;
746 }
747
755 for (int row = 0; row < rows(); ++row) {
756 for (int col = 0; col < rhs.cols(); ++col) {
757 (*this)[row, col] *= rhs;
758 }
759 }
760
761 return *this;
762 }
763
771 template <EigenMatrixLike LHS>
773 const Variable<Scalar>& rhs) {
774 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
775
776 for (int row = 0; row < result.rows(); ++row) {
777 for (int col = 0; col < result.cols(); ++col) {
778 result[row, col] = lhs(row, col) / rhs;
779 }
780 }
781
782 return result;
783 }
784
792 template <SleipnirMatrixLike<Scalar> LHS, ScalarLike RHS>
794 friend VariableMatrix<Scalar> operator/(const LHS& lhs, const RHS& rhs) {
795 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
796
797 for (int row = 0; row < result.rows(); ++row) {
798 for (int col = 0; col < result.cols(); ++col) {
799 result[row, col] = lhs[row, col] / rhs;
800 }
801 }
802
803 return result;
804 }
805
813 template <SleipnirMatrixLike<Scalar> LHS>
815 const Variable<Scalar>& rhs) {
816 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
817
818 for (int row = 0; row < result.rows(); ++row) {
819 for (int col = 0; col < result.cols(); ++col) {
820 result[row, col] = lhs[row, col] / rhs;
821 }
822 }
823
824 return result;
825 }
826
834 for (int row = 0; row < rows(); ++row) {
835 for (int col = 0; col < cols(); ++col) {
836 (*this)[row, col] /= rhs;
837 }
838 }
839
840 return *this;
841 }
842
850 template <EigenMatrixLike LHS, SleipnirMatrixLike<Scalar> RHS>
852 slp_assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols());
853
854 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
855
856 for (int row = 0; row < result.rows(); ++row) {
857 for (int col = 0; col < result.cols(); ++col) {
858 result[row, col] = lhs(row, col) + rhs[row, col];
859 }
860 }
861
862 return result;
863 }
864
872 template <SleipnirMatrixLike<Scalar> LHS, EigenMatrixLike RHS>
874 slp_assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols());
875
876 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
877
878 for (int row = 0; row < result.rows(); ++row) {
879 for (int col = 0; col < result.cols(); ++col) {
880 result[row, col] = lhs[row, col] + rhs(row, col);
881 }
882 }
883
884 return result;
885 }
886
894 template <SleipnirMatrixLike<Scalar> LHS, SleipnirMatrixLike<Scalar> RHS>
896 slp_assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols());
897
898 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
899
900 for (int row = 0; row < result.rows(); ++row) {
901 for (int col = 0; col < result.cols(); ++col) {
902 result[row, col] = lhs[row, col] + rhs[row, col];
903 }
904 }
905
906 return result;
907 }
908
916 slp_assert(rows() == rhs.rows() && cols() == rhs.cols());
917
918 for (int row = 0; row < rows(); ++row) {
919 for (int col = 0; col < cols(); ++col) {
920 if constexpr (SleipnirMatrixLike<decltype(rhs), Scalar>) {
921 (*this)[row, col] += rhs[row, col];
922 } else {
923 (*this)[row, col] += rhs(row, col);
924 }
925 }
926 }
927
928 return *this;
929 }
930
938 slp_assert(rows() == 1 && cols() == 1);
939
940 for (int row = 0; row < rows(); ++row) {
941 for (int col = 0; col < cols(); ++col) {
942 (*this)[row, col] += rhs;
943 }
944 }
945
946 return *this;
947 }
948
956 template <EigenMatrixLike LHS, SleipnirMatrixLike<Scalar> RHS>
958 slp_assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols());
959
960 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
961
962 for (int row = 0; row < result.rows(); ++row) {
963 for (int col = 0; col < result.cols(); ++col) {
964 result[row, col] = lhs(row, col) - rhs[row, col];
965 }
966 }
967
968 return result;
969 }
970
978 template <SleipnirMatrixLike<Scalar> LHS, EigenMatrixLike RHS>
980 slp_assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols());
981
982 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
983
984 for (int row = 0; row < result.rows(); ++row) {
985 for (int col = 0; col < result.cols(); ++col) {
986 result[row, col] = lhs[row, col] - rhs(row, col);
987 }
988 }
989
990 return result;
991 }
992
1000 template <SleipnirMatrixLike<Scalar> LHS, SleipnirMatrixLike<Scalar> RHS>
1002 slp_assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols());
1003
1004 VariableMatrix<Scalar> result(detail::empty, lhs.rows(), lhs.cols());
1005
1006 for (int row = 0; row < result.rows(); ++row) {
1007 for (int col = 0; col < result.cols(); ++col) {
1008 result[row, col] = lhs[row, col] - rhs[row, col];
1009 }
1010 }
1011
1012 return result;
1013 }
1014
1022 slp_assert(rows() == rhs.rows() && cols() == rhs.cols());
1023
1024 for (int row = 0; row < rows(); ++row) {
1025 for (int col = 0; col < cols(); ++col) {
1026 if constexpr (SleipnirMatrixLike<decltype(rhs), Scalar>) {
1027 (*this)[row, col] -= rhs[row, col];
1028 } else {
1029 (*this)[row, col] -= rhs(row, col);
1030 }
1031 }
1032 }
1033
1034 return *this;
1035 }
1036
1044 slp_assert(rows() == 1 && cols() == 1);
1045
1046 for (int row = 0; row < rows(); ++row) {
1047 for (int col = 0; col < cols(); ++col) {
1048 (*this)[row, col] -= rhs;
1049 }
1050 }
1051
1052 return *this;
1053 }
1054
1061 const SleipnirMatrixLike<Scalar> auto& lhs) {
1062 VariableMatrix<Scalar> result{detail::empty, lhs.rows(), lhs.cols()};
1063
1064 for (int row = 0; row < result.rows(); ++row) {
1065 for (int col = 0; col < result.cols(); ++col) {
1066 result[row, col] = -lhs[row, col];
1067 }
1068 }
1069
1070 return result;
1071 }
1072
1076 // NOLINTNEXTLINE (google-explicit-constructor)
1077 operator Variable<Scalar>() const {
1078 slp_assert(rows() == 1 && cols() == 1);
1079 return (*this)[0, 0];
1080 }
1081
1088 VariableMatrix<Scalar> result{detail::empty, cols(), rows()};
1089
1090 for (int row = 0; row < rows(); ++row) {
1091 for (int col = 0; col < cols(); ++col) {
1092 result[col, row] = (*this)[row, col];
1093 }
1094 }
1095
1096 return result;
1097 }
1098
1104 int rows() const { return m_rows; }
1105
1111 int cols() const { return m_cols; }
1112
1120 Scalar value(int row, int col) { return (*this)[row, col].value(); }
1121
1128 Scalar value(int index) { return (*this)[index].value(); }
1129
1135 Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> value() {
1136 Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> result{rows(),
1137 cols()};
1138
1139 for (int row = 0; row < rows(); ++row) {
1140 for (int col = 0; col < cols(); ++col) {
1141 result(row, col) = value(row, col);
1142 }
1143 }
1144
1145 return result;
1146 }
1147
1156 const {
1157 VariableMatrix<Scalar> result{detail::empty, rows(), cols()};
1158
1159 for (int row = 0; row < rows(); ++row) {
1160 for (int col = 0; col < cols(); ++col) {
1161 result[row, col] = unary_op((*this)[row, col]);
1162 }
1163 }
1164
1165 return result;
1166 }
1167
1168#ifndef DOXYGEN_SHOULD_SKIP_THIS
1169
1170 class iterator {
1171 public:
1172 using iterator_category = std::bidirectional_iterator_tag;
1173 using value_type = Variable<Scalar>;
1174 using difference_type = std::ptrdiff_t;
1175 using pointer = Variable<Scalar>*;
1176 using reference = Variable<Scalar>&;
1177
1178 constexpr iterator() noexcept = default;
1179
1182 : m_it{it} {}
1183
1184 constexpr iterator& operator++() noexcept {
1185 ++m_it;
1186 return *this;
1187 }
1188
1189 constexpr iterator operator++(int) noexcept {
1190 iterator retval = *this;
1191 ++(*this);
1192 return retval;
1193 }
1194
1195 constexpr iterator& operator--() noexcept {
1196 --m_it;
1197 return *this;
1198 }
1199
1200 constexpr iterator operator--(int) noexcept {
1201 iterator retval = *this;
1202 --(*this);
1203 return retval;
1204 }
1205
1206 constexpr bool operator==(const iterator&) const noexcept = default;
1207
1208 constexpr reference operator*() const noexcept { return *m_it; }
1209
1210 private:
1211 gch::small_vector<Variable<Scalar>>::iterator m_it;
1212 };
1213
1214 class const_iterator {
1215 public:
1216 using iterator_category = std::bidirectional_iterator_tag;
1217 using value_type = Variable<Scalar>;
1218 using difference_type = std::ptrdiff_t;
1219 using pointer = Variable<Scalar>*;
1220 using const_reference = const Variable<Scalar>&;
1221
1222 constexpr const_iterator() noexcept = default;
1223
1224 explicit constexpr const_iterator(
1225 gch::small_vector<Variable<Scalar>>::const_iterator it) noexcept
1226 : m_it{it} {}
1227
1228 constexpr const_iterator& operator++() noexcept {
1229 ++m_it;
1230 return *this;
1231 }
1232
1233 constexpr const_iterator operator++(int) noexcept {
1234 const_iterator retval = *this;
1235 ++(*this);
1236 return retval;
1237 }
1238
1239 constexpr const_iterator& operator--() noexcept {
1240 --m_it;
1241 return *this;
1242 }
1243
1244 constexpr const_iterator operator--(int) noexcept {
1245 const_iterator retval = *this;
1246 --(*this);
1247 return retval;
1248 }
1249
1250 constexpr bool operator==(const const_iterator&) const noexcept = default;
1251
1252 constexpr const_reference operator*() const noexcept { return *m_it; }
1253
1254 private:
1255 gch::small_vector<Variable<Scalar>>::const_iterator m_it;
1256 };
1257
1258 using reverse_iterator = std::reverse_iterator<iterator>;
1259 using const_reverse_iterator = std::reverse_iterator<const_iterator>;
1260
1261#endif // DOXYGEN_SHOULD_SKIP_THIS
1262
1268 iterator begin() { return iterator{m_storage.begin()}; }
1269
1275 iterator end() { return iterator{m_storage.end()}; }
1276
1282 const_iterator begin() const { return const_iterator{m_storage.begin()}; }
1283
1289 const_iterator end() const { return const_iterator{m_storage.end()}; }
1290
1296 const_iterator cbegin() const { return const_iterator{m_storage.cbegin()}; }
1297
1303 const_iterator cend() const { return const_iterator{m_storage.cend()}; }
1304
1311
1318
1327
1336
1345
1354
1360 size_t size() const { return m_storage.size(); }
1361
1370 VariableMatrix<Scalar> result{detail::empty, rows, cols};
1371
1372 for (auto& elem : result) {
1373 elem = Scalar(0);
1374 }
1375
1376 return result;
1377 }
1378
1387 VariableMatrix<Scalar> result{detail::empty, rows, cols};
1388
1389 for (auto& elem : result) {
1390 elem = Scalar(1);
1391 }
1392
1393 return result;
1394 }
1395
1396 private:
1397 gch::small_vector<Variable<Scalar>> m_storage;
1398 int m_rows = 0;
1399 int m_cols = 0;
1400};
1401
1402template <typename Derived>
1403VariableMatrix(const Eigen::MatrixBase<Derived>&)
1404 -> VariableMatrix<typename Derived::Scalar>;
1405
1406template <typename Derived>
1407VariableMatrix(const Eigen::DiagonalBase<Derived>&)
1408 -> VariableMatrix<typename Derived::Scalar>;
1409
1418template <typename Scalar>
1419VariableMatrix<Scalar> cwise_reduce(
1420 const VariableMatrix<Scalar>& lhs, const VariableMatrix<Scalar>& rhs,
1421 function_ref<Variable<Scalar>(const Variable<Scalar>& x,
1422 const Variable<Scalar>& y)>
1423 binary_op) {
1424 slp_assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols());
1425
1426 VariableMatrix<Scalar> result{detail::empty, lhs.rows(), lhs.cols()};
1427
1428 for (int row = 0; row < lhs.rows(); ++row) {
1429 for (int col = 0; col < lhs.cols(); ++col) {
1430 result[row, col] = binary_op(lhs[row, col], rhs[row, col]);
1431 }
1432 }
1433
1434 return result;
1435}
1436
1448template <typename Scalar>
1449VariableMatrix<Scalar> block(
1450 std::initializer_list<std::initializer_list<VariableMatrix<Scalar>>> list) {
1451 // Get row and column counts for destination matrix
1452 int rows = 0;
1453 int cols = -1;
1454 for (const auto& row : list) {
1455 if (row.size() > 0) {
1456 rows += row.begin()->rows();
1457 }
1458
1459 // Get number of columns in this row
1460 int latest_cols = 0;
1461 for (const auto& elem : row) {
1462 // Assert the first and latest row have the same height
1463 slp_assert(row.begin()->rows() == elem.rows());
1464
1465 latest_cols += elem.cols();
1466 }
1467
1468 // If this is the first row, record the column count. Otherwise, assert the
1469 // first and latest column counts are the same.
1470 if (cols == -1) {
1471 cols = latest_cols;
1472 } else {
1473 slp_assert(cols == latest_cols);
1474 }
1475 }
1476
1477 VariableMatrix<Scalar> result{detail::empty, rows, cols};
1478
1479 int row_offset = 0;
1480 for (const auto& row : list) {
1481 int col_offset = 0;
1482 for (const auto& elem : row) {
1483 result.block(row_offset, col_offset, elem.rows(), elem.cols()) = elem;
1484 col_offset += elem.cols();
1485 }
1486 if (row.size() > 0) {
1487 row_offset += row.begin()->rows();
1488 }
1489 }
1490
1491 return result;
1492}
1493
1507template <typename Scalar>
1508VariableMatrix<Scalar> block(
1509 const std::vector<std::vector<VariableMatrix<Scalar>>>& list) {
1510 // Get row and column counts for destination matrix
1511 int rows = 0;
1512 int cols = -1;
1513 for (const auto& row : list) {
1514 if (row.size() > 0) {
1515 rows += row.begin()->rows();
1516 }
1517
1518 // Get number of columns in this row
1519 int latest_cols = 0;
1520 for (const auto& elem : row) {
1521 // Assert the first and latest row have the same height
1522 slp_assert(row.begin()->rows() == elem.rows());
1523
1524 latest_cols += elem.cols();
1525 }
1526
1527 // If this is the first row, record the column count. Otherwise, assert the
1528 // first and latest column counts are the same.
1529 if (cols == -1) {
1530 cols = latest_cols;
1531 } else {
1532 slp_assert(cols == latest_cols);
1533 }
1534 }
1535
1536 VariableMatrix<Scalar> result{detail::empty, rows, cols};
1537
1538 int row_offset = 0;
1539 for (const auto& row : list) {
1540 int col_offset = 0;
1541 for (const auto& elem : row) {
1542 result.block(row_offset, col_offset, elem.rows(), elem.cols()) = elem;
1543 col_offset += elem.cols();
1544 }
1545 if (row.size() > 0) {
1546 row_offset += row.begin()->rows();
1547 }
1548 }
1549
1550 return result;
1551}
1552
1561template <typename Scalar>
1562VariableMatrix<Scalar> solve(const VariableMatrix<Scalar>& A,
1563 const VariableMatrix<Scalar>& B) {
1564 // m x n * n x p = m x p
1565 slp_assert(A.rows() == B.rows());
1566
1567 if (A.rows() == 1 && A.cols() == 1) {
1568 // Compute optimal inverse instead of using Eigen's general solver
1569 return B[0, 0] / A[0, 0];
1570 } else if (A.rows() == 2 && A.cols() == 2) {
1571 // Compute optimal inverse instead of using Eigen's general solver
1572 //
1573 // [a b]⁻¹ ___1___ [ d −b]
1574 // [c d] = ad − bc [−c a]
1575
1576 const auto& a = A[0, 0];
1577 const auto& b = A[0, 1];
1578 const auto& c = A[1, 0];
1579 const auto& d = A[1, 1];
1580
1581 VariableMatrix adj_A{{d, -b}, {-c, a}};
1582 auto det_A = a * d - b * c;
1583 return adj_A / det_A * B;
1584 } else if (A.rows() == 3 && A.cols() == 3) {
1585 // Compute optimal inverse instead of using Eigen's general solver
1586 //
1587 // [a b c]⁻¹
1588 // [d e f]
1589 // [g h i]
1590 // 1 [ei − fh ch − bi bf − ce]
1591 // = ------------------------------------ [fg − di ai − cg cd − af]
1592 // a(ei − fh) + b(fg − di) + c(dh − eg) [dh − eg bg − ah ae − bd]
1593 //
1594 // https://www.wolframalpha.com/input?i=inverse+%7B%7Ba%2C+b%2C+c%7D%2C+%7Bd%2C+e%2C+f%7D%2C+%7Bg%2C+h%2C+i%7D%7D
1595
1596 const auto& a = A[0, 0];
1597 const auto& b = A[0, 1];
1598 const auto& c = A[0, 2];
1599 const auto& d = A[1, 0];
1600 const auto& e = A[1, 1];
1601 const auto& f = A[1, 2];
1602 const auto& g = A[2, 0];
1603 const auto& h = A[2, 1];
1604 const auto& i = A[2, 2];
1605
1606 auto ae = a * e;
1607 auto af = a * f;
1608 auto ah = a * h;
1609 auto ai = a * i;
1610 auto bd = b * d;
1611 auto bf = b * f;
1612 auto bg = b * g;
1613 auto bi = b * i;
1614 auto cd = c * d;
1615 auto ce = c * e;
1616 auto cg = c * g;
1617 auto ch = c * h;
1618 auto dh = d * h;
1619 auto di = d * i;
1620 auto eg = e * g;
1621 auto ei = e * i;
1622 auto fg = f * g;
1623 auto fh = f * h;
1624
1625 auto adj_A00 = ei - fh;
1626 auto adj_A10 = fg - di;
1627 auto adj_A20 = dh - eg;
1628
1629 VariableMatrix adj_A{{adj_A00, ch - bi, bf - ce},
1630 {adj_A10, ai - cg, cd - af},
1631 {adj_A20, bg - ah, ae - bd}};
1632 auto det_A = a * adj_A00 + b * adj_A10 + c * adj_A20;
1633 return adj_A / det_A * B;
1634 } else if (A.rows() == 4 && A.cols() == 4) {
1635 // Compute optimal inverse instead of using Eigen's general solver
1636 //
1637 // [a b c d]⁻¹
1638 // [e f g h]
1639 // [i j k l]
1640 // [m n o p]
1641 //
1642 // https://www.wolframalpha.com/input?i=inverse+%7B%7Ba%2C+b%2C+c%2C+d%7D%2C+%7Be%2C+f%2C+g%2C+h%7D%2C+%7Bi%2C+j%2C+k%2C+l%7D%2C+%7Bm%2C+n%2C+o%2C+p%7D%7D
1643
1644 const auto& a = A[0, 0];
1645 const auto& b = A[0, 1];
1646 const auto& c = A[0, 2];
1647 const auto& d = A[0, 3];
1648 const auto& e = A[1, 0];
1649 const auto& f = A[1, 1];
1650 const auto& g = A[1, 2];
1651 const auto& h = A[1, 3];
1652 const auto& i = A[2, 0];
1653 const auto& j = A[2, 1];
1654 const auto& k = A[2, 2];
1655 const auto& l = A[2, 3];
1656 const auto& m = A[3, 0];
1657 const auto& n = A[3, 1];
1658 const auto& o = A[3, 2];
1659 const auto& p = A[3, 3];
1660
1661 auto afk = a * f * k;
1662 auto afl = a * f * l;
1663 auto afo = a * f * o;
1664 auto afp = a * f * p;
1665 auto agj = a * g * j;
1666 auto agl = a * g * l;
1667 auto agn = a * g * n;
1668 auto agp = a * g * p;
1669 auto ahj = a * h * j;
1670 auto ahk = a * h * k;
1671 auto ahn = a * h * n;
1672 auto aho = a * h * o;
1673 auto ajo = a * j * o;
1674 auto ajp = a * j * p;
1675 auto akn = a * k * n;
1676 auto akp = a * k * p;
1677 auto aln = a * l * n;
1678 auto alo = a * l * o;
1679 auto bek = b * e * k;
1680 auto bel = b * e * l;
1681 auto beo = b * e * o;
1682 auto bep = b * e * p;
1683 auto bgi = b * g * i;
1684 auto bgl = b * g * l;
1685 auto bgm = b * g * m;
1686 auto bgp = b * g * p;
1687 auto bhi = b * h * i;
1688 auto bhk = b * h * k;
1689 auto bhm = b * h * m;
1690 auto bho = b * h * o;
1691 auto bio = b * i * o;
1692 auto bip = b * i * p;
1693 auto bjp = b * j * p;
1694 auto bkm = b * k * m;
1695 auto bkp = b * k * p;
1696 auto blm = b * l * m;
1697 auto blo = b * l * o;
1698 auto cej = c * e * j;
1699 auto cel = c * e * l;
1700 auto cen = c * e * n;
1701 auto cep = c * e * p;
1702 auto cfi = c * f * i;
1703 auto cfl = c * f * l;
1704 auto cfm = c * f * m;
1705 auto cfp = c * f * p;
1706 auto chi = c * h * i;
1707 auto chj = c * h * j;
1708 auto chm = c * h * m;
1709 auto chn = c * h * n;
1710 auto cin = c * i * n;
1711 auto cip = c * i * p;
1712 auto cjm = c * j * m;
1713 auto cjp = c * j * p;
1714 auto clm = c * l * m;
1715 auto cln = c * l * n;
1716 auto dej = d * e * j;
1717 auto dek = d * e * k;
1718 auto den = d * e * n;
1719 auto deo = d * e * o;
1720 auto dfi = d * f * i;
1721 auto dfk = d * f * k;
1722 auto dfm = d * f * m;
1723 auto dfo = d * f * o;
1724 auto dgi = d * g * i;
1725 auto dgj = d * g * j;
1726 auto dgm = d * g * m;
1727 auto dgn = d * g * n;
1728 auto din = d * i * n;
1729 auto dio = d * i * o;
1730 auto djm = d * j * m;
1731 auto djo = d * j * o;
1732 auto dkm = d * k * m;
1733 auto dkn = d * k * n;
1734 auto ejo = e * j * o;
1735 auto ejp = e * j * p;
1736 auto ekn = e * k * n;
1737 auto ekp = e * k * p;
1738 auto eln = e * l * n;
1739 auto elo = e * l * o;
1740 auto fio = f * i * o;
1741 auto fip = f * i * p;
1742 auto fkm = f * k * m;
1743 auto fkp = f * k * p;
1744 auto flm = f * l * m;
1745 auto flo = f * l * o;
1746 auto gin = g * i * n;
1747 auto gip = g * i * p;
1748 auto gjm = g * j * m;
1749 auto gjp = g * j * p;
1750 auto glm = g * l * m;
1751 auto gln = g * l * n;
1752 auto hin = h * i * n;
1753 auto hio = h * i * o;
1754 auto hjm = h * j * m;
1755 auto hjo = h * j * o;
1756 auto hkm = h * k * m;
1757 auto hkn = h * k * n;
1758
1759 auto adj_A00 = fkp - flo - gjp + gln + hjo - hkn;
1760 auto adj_A01 = -bkp + blo + cjp - cln - djo + dkn;
1761 auto adj_A02 = bgp - bho - cfp + chn + dfo - dgn;
1762 auto adj_A03 = -bgl + bhk + cfl - chj - dfk + dgj;
1763 auto adj_A10 = -ekp + elo + gip - glm - hio + hkm;
1764 auto adj_A11 = akp - alo - cip + clm + dio - dkm;
1765 auto adj_A12 = -agp + aho + cep - chm - deo + dgm;
1766 auto adj_A13 = agl - ahk - cel + chi + dek - dgi;
1767 auto adj_A20 = ejp - eln - fip + flm + hin - hjm;
1768 auto adj_A21 = -ajp + aln + bip - blm - din + djm;
1769 auto adj_A22 = afp - ahn - bep + bhm + den - dfm;
1770 auto adj_A23 = -afl + ahj + bel - bhi - dej + dfi;
1771 auto adj_A30 = -ejo + ekn + fio - fkm - gin + gjm;
1772 // NOLINTNEXTLINE(build/include_what_you_use)
1773 auto adj_A31 = ajo - akn - bio + bkm + cin - cjm;
1774 auto adj_A32 = -afo + agn + beo - bgm - cen + cfm;
1775 auto adj_A33 = afk - agj - bek + bgi + cej - cfi;
1776
1777 VariableMatrix adj_A{{adj_A00, adj_A01, adj_A02, adj_A03},
1778 {adj_A10, adj_A11, adj_A12, adj_A13},
1779 {adj_A20, adj_A21, adj_A22, adj_A23},
1780 {adj_A30, adj_A31, adj_A32, adj_A33}};
1781 auto det_A = a * adj_A00 + b * adj_A10 + c * adj_A20 + d * adj_A30;
1782 return adj_A / det_A * B;
1783 } else {
1784 using MatrixXv =
1785 Eigen::Matrix<Variable<Scalar>, Eigen::Dynamic, Eigen::Dynamic>;
1786
1787 MatrixXv eigen_A{A.rows(), A.cols()};
1788 for (int row = 0; row < A.rows(); ++row) {
1789 for (int col = 0; col < A.cols(); ++col) {
1790 eigen_A(row, col) = A[row, col];
1791 }
1792 }
1793
1794 MatrixXv eigen_B{B.rows(), B.cols()};
1795 for (int row = 0; row < B.rows(); ++row) {
1796 for (int col = 0; col < B.cols(); ++col) {
1797 eigen_B(row, col) = B[row, col];
1798 }
1799 }
1800
1801 MatrixXv eigen_X = eigen_A.householderQr().solve(eigen_B);
1802
1803 VariableMatrix<Scalar> X{detail::empty, A.cols(), B.cols()};
1804 for (int row = 0; row < X.rows(); ++row) {
1805 for (int col = 0; col < X.cols(); ++col) {
1806 X[row, col] = eigen_X(row, col);
1807 }
1808 }
1809
1810 return X;
1811 }
1812}
1813
1814extern template SLEIPNIR_DLLEXPORT VariableMatrix<double> solve(
1815 const VariableMatrix<double>& A, const VariableMatrix<double>& B);
1816
1817} // namespace slp
Definition intrusive_shared_ptr.hpp:29
Definition sleipnir_base.hpp:11
Definition slice.hpp:31
Definition variable_block.hpp:26
Definition variable_matrix.hpp:35
const_reverse_iterator crbegin() const
Definition variable_matrix.hpp:1342
VariableMatrix(std::initializer_list< std::initializer_list< Variable< Scalar > > > list)
Definition variable_matrix.hpp:87
const_reverse_iterator crend() const
Definition variable_matrix.hpp:1351
Scalar_ Scalar
Definition variable_matrix.hpp:40
iterator end()
Definition variable_matrix.hpp:1275
const Variable< Scalar > & operator[](int row, int col) const
Definition variable_matrix.hpp:361
VariableMatrix(Variable< Scalar > &&variable)
Definition variable_matrix.hpp:220
const VariableBlock< const VariableMatrix > block(int row_offset, int col_offset, int block_rows, int block_cols) const
Definition variable_matrix.hpp:416
size_t size() const
Definition variable_matrix.hpp:1360
static VariableMatrix< Scalar > zero(int rows, int cols)
Definition variable_matrix.hpp:1369
VariableBlock< VariableMatrix > block(int row_offset, int col_offset, int block_rows, int block_cols)
Definition variable_matrix.hpp:398
VariableMatrix(const std::vector< std::vector< Scalar > > &list)
Definition variable_matrix.hpp:116
VariableMatrix & operator=(ScalarLike auto value)
Definition variable_matrix.hpp:316
VariableMatrix & operator-=(const MatrixLike auto &rhs)
Definition variable_matrix.hpp:1021
Variable< Scalar > & operator[](int index)
Definition variable_matrix.hpp:373
VariableMatrix & operator=(const Eigen::MatrixBase< Derived > &values)
Definition variable_matrix.hpp:296
friend VariableMatrix< Scalar > operator-(const LHS &lhs, const RHS &rhs)
Definition variable_matrix.hpp:957
VariableMatrix & operator-=(const ScalarLike auto &rhs)
Definition variable_matrix.hpp:1043
friend VariableMatrix< Scalar > operator*(const Variable< Scalar > &lhs, const RHS &rhs)
Definition variable_matrix.hpp:690
VariableMatrix & operator*=(const ScalarLike auto &rhs)
Definition variable_matrix.hpp:754
Scalar value(int index)
Definition variable_matrix.hpp:1128
VariableMatrix(const VariableBlock< VariableMatrix > &values)
Definition variable_matrix.hpp:230
VariableMatrix(int rows, int cols)
Definition variable_matrix.hpp:61
const_iterator cbegin() const
Definition variable_matrix.hpp:1296
const_iterator begin() const
Definition variable_matrix.hpp:1282
const_iterator cend() const
Definition variable_matrix.hpp:1303
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > value()
Definition variable_matrix.hpp:1135
const VariableBlock< const VariableMatrix > operator[](Slice row_slice, Slice col_slice) const
Definition variable_matrix.hpp:448
friend VariableMatrix< Scalar > operator/(const LHS &lhs, const Variable< Scalar > &rhs)
Definition variable_matrix.hpp:772
VariableMatrix(int rows)
Definition variable_matrix.hpp:53
VariableMatrix(std::span< const Variable< Scalar > > values, int rows, int cols)
Definition variable_matrix.hpp:278
VariableBlock< VariableMatrix > col(int col)
Definition variable_matrix.hpp:553
VariableMatrix(const Eigen::MatrixBase< Derived > &values)
Definition variable_matrix.hpp:171
const VariableBlock< const VariableMatrix > operator[](Slice row_slice, int row_slice_length, Slice col_slice, int col_slice_length) const
Definition variable_matrix.hpp:489
const VariableBlock< const VariableMatrix > col(int col) const
Definition variable_matrix.hpp:564
VariableMatrix & operator+=(const ScalarLike auto &rhs)
Definition variable_matrix.hpp:937
friend VariableMatrix< Scalar > operator*(const LHS &lhs, const Variable< Scalar > &rhs)
Definition variable_matrix.hpp:651
VariableBlock< VariableMatrix > segment(int offset, int length)
Definition variable_matrix.hpp:503
Scalar value(int row, int col)
Definition variable_matrix.hpp:1120
const VariableBlock< const VariableMatrix > row(int row) const
Definition variable_matrix.hpp:542
VariableMatrix()=default
static VariableMatrix< Scalar > ones(int rows, int cols)
Definition variable_matrix.hpp:1386
VariableMatrix & operator+=(const MatrixLike auto &rhs)
Definition variable_matrix.hpp:915
VariableMatrix(std::span< const Variable< Scalar > > values)
Definition variable_matrix.hpp:261
friend VariableMatrix< Scalar > operator-(const SleipnirMatrixLike< Scalar > auto &lhs)
Definition variable_matrix.hpp:1060
VariableMatrix< Scalar > T() const
Definition variable_matrix.hpp:1087
VariableMatrix< Scalar > cwise_transform(function_ref< Variable< Scalar >(const Variable< Scalar > &x)> unary_op) const
Definition variable_matrix.hpp:1154
friend VariableMatrix< Scalar > operator*(const LHS &lhs, const RHS &rhs)
Definition variable_matrix.hpp:576
const_reverse_iterator rbegin() const
Definition variable_matrix.hpp:1324
VariableMatrix(const Variable< Scalar > &variable)
Definition variable_matrix.hpp:210
reverse_iterator rend()
Definition variable_matrix.hpp:1317
VariableMatrix(const VariableBlock< const VariableMatrix > &values)
Definition variable_matrix.hpp:246
int rows() const
Definition variable_matrix.hpp:1104
VariableMatrix(const Eigen::DiagonalBase< Derived > &values)
Definition variable_matrix.hpp:189
VariableBlock< VariableMatrix > operator[](Slice row_slice, Slice col_slice)
Definition variable_matrix.hpp:434
const_reverse_iterator rend() const
Definition variable_matrix.hpp:1333
VariableMatrix & operator*=(const MatrixLike auto &rhs)
Definition variable_matrix.hpp:728
VariableMatrix(detail::empty_t, int rows, int cols)
Definition variable_matrix.hpp:74
const_iterator end() const
Definition variable_matrix.hpp:1289
void set_value(const Eigen::MatrixBase< Derived > &values)
Definition variable_matrix.hpp:331
Variable< Scalar > & operator[](int row, int col)
Definition variable_matrix.hpp:348
iterator begin()
Definition variable_matrix.hpp:1268
const VariableBlock< const VariableMatrix > segment(int offset, int length) const
Definition variable_matrix.hpp:517
reverse_iterator rbegin()
Definition variable_matrix.hpp:1310
VariableBlock< VariableMatrix > row(int row)
Definition variable_matrix.hpp:531
friend VariableMatrix< Scalar > operator+(const LHS &lhs, const RHS &rhs)
Definition variable_matrix.hpp:851
VariableMatrix(const std::vector< std::vector< Variable< Scalar > > > &list)
Definition variable_matrix.hpp:144
VariableBlock< VariableMatrix > operator[](Slice row_slice, int row_slice_length, Slice col_slice, int col_slice_length)
Definition variable_matrix.hpp:469
VariableMatrix & operator/=(const ScalarLike auto &rhs)
Definition variable_matrix.hpp:833
const Variable< Scalar > & operator[](int index) const
Definition variable_matrix.hpp:384
int cols() const
Definition variable_matrix.hpp:1111
Definition variable.hpp:49
Definition function_ref.hpp:13
Definition concepts.hpp:18
Definition concepts.hpp:24
Definition concepts.hpp:33
Definition concepts.hpp:38
Definition empty.hpp:10