proxsuite 0.7.2
The Advanced Proximal Optimization Toolbox
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Loading...
Searching...
No Matches
views.hpp
Go to the documentation of this file.
1//
2// Copyright (c) 2022 INRIA
3//
7#ifndef PROXSUITE_PROXQP_DENSE_VIEWS_HPP
8#define PROXSUITE_PROXQP_DENSE_VIEWS_HPP
9
12#include <cstring>
13#include <type_traits>
14#include <Eigen/Core>
15
16#define LDLT_CONCEPT(...) \
17 VEG_CONCEPT_MACRO(::proxsuite::proxqp::concepts, __VA_ARGS__)
18#define LDLT_CHECK_CONCEPT(...) \
19 VEG_CHECK_CONCEPT_MACRO(::proxqp::concepts, __VA_ARGS__)
20
21namespace proxsuite {
22namespace proxqp {
23
24using usize = decltype(sizeof(0));
25namespace detail {
26template<typename Fn>
27struct FnInfo;
28template<typename Ret_, typename... Args>
29struct FnInfo<auto(Args...)->Ret_>
30{
31 template<usize I>
32 using Arg = proxsuite::linalg::veg::ith<I, Args...>;
33 using Ret = Ret_;
34};
35} // namespace detail
36
37#define LDLT_IMPL_GET_PARAM(Fn, Idx) \
38 typename ::proxsuite::proxqp::detail::FnInfo< \
39 decltype Fn /* NOLINT */>::template Arg<(Idx)>,
40
41#define LDLT_IMPL_GET_PARAMS_0(NParams, ...) \
42 __VEG_PP_TUPLE_FOR_EACH(LDLT_IMPL_GET_PARAM, \
43 (__VA_ARGS__), \
44 __VEG_PP_MAKE_TUPLE(__VEG_IMPL_PP_DEC(NParams)))
45
46#define LDLT_IMPL_GET_PARAMS_1(NParams, ...)
47
48#define LDLT_IMPL_GET_PARAMS(NParams, ...) \
49 __VEG_PP_CAT2(LDLT_IMPL_GET_PARAMS_, __VEG_IMPL_PP_IS_1(NParams)) \
50 (NParams, __VA_ARGS__)
51
52#define LDLT_EXPLICIT_TPL_DEF(NParams, ...) \
53 template auto __VA_ARGS__( \
54 LDLT_IMPL_GET_PARAMS(NParams, __VA_ARGS__) \
55 typename ::proxsuite::proxqp::detail::FnInfo< \
56 decltype(__VA_ARGS__)>::template Arg<(NParams) - 1>) -> \
57 typename ::proxsuite::proxqp::detail::FnInfo<decltype(__VA_ARGS__)>::Ret
58#define LDLT_EXPLICIT_TPL_DECL(NParams, ...) \
59 extern LDLT_EXPLICIT_TPL_DEF(NParams, __VA_ARGS__)
60
67using f32 = float;
68using f64 = double;
69namespace detail {
70
71struct NoCopy
72{
73 NoCopy() = default;
74 ~NoCopy() = default;
75
76 NoCopy(NoCopy const&) = delete;
77 NoCopy(NoCopy&&) = delete;
78 auto operator=(NoCopy const&) -> NoCopy& = delete;
79 auto operator=(NoCopy&&) -> NoCopy& = delete;
80};
81
82template<typename Fn>
83struct Defer /* NOLINT */
84{
85 Fn fn;
87
88 VEG_INLINE ~Defer() noexcept(noexcept(VEG_FWD(fn)())) { VEG_FWD(fn)(); }
89};
90
91namespace nb {
92struct defer
93{
94 template<typename Fn>
95 VEG_INLINE constexpr auto operator()(Fn fn) const -> Defer<Fn>
96 {
97 return { VEG_FWD(fn), {} };
98 }
99};
100struct max2
101{
102 template<typename T>
103 VEG_INLINE constexpr auto operator()(T const& a, T const& b) const -> T const&
104 {
105 return a > b ? a : b;
106 }
107};
108struct min2
109{
110 template<typename T>
111 VEG_INLINE constexpr auto operator()(T a, T b) const -> T
112 {
113 return (a < b) ? a : b;
114 }
115};
116
117} // namespace nb
118
119template<typename T>
120constexpr auto
121min_list_impl(T init, T const* arr, usize n) noexcept -> T
122{
123 return (n == 0)
124 ? init
125 : nb::min2{}(init, detail::min_list_impl(*arr, arr + 1, n - 1));
126}
127template<typename T, usize N>
128constexpr auto
129cx_min_list(T const (&arr)[N]) noexcept -> T
130{
131 return detail::min_list_impl( //
132 arr[0],
133 arr + 1,
134 N - 1);
135}
136
137namespace nb {
139{
140 template<typename T>
141 VEG_INLINE auto operator()(std::initializer_list<T> list) const -> T
142 {
143 T const* data = list.begin();
144 isize len = isize(list.size());
145
146 T current_max = data[0];
147 for (isize i = 1; i < len; ++i) {
148 if (data[i] > current_max) {
149 current_max = data[i];
150 }
151 }
152 return current_max;
153 }
154};
155} // namespace nb
159VEG_NIEBLOID(max_list);
160
161template<typename T, bool = std::is_floating_point<T>::value>
163{
164 static void fn(T* dest, usize n)
165 {
166 for (usize i = 0; i < n; ++i) {
167 *dest = 0;
168 }
169 }
170};
171
172template<typename T>
173struct SetZeroImpl<T, true>
174{
175 static void fn(T* dest, usize n)
176 {
177 // TODO: assert bit representation is zero
178 std::memset(dest, 0, n * sizeof(T));
179 }
180};
181
182template<typename T>
183void
184set_zero(T* dest, usize n)
185{
186 SetZeroImpl<T>::fn(dest, n);
187}
188
189constexpr auto
190round_up(isize n, isize k) noexcept -> isize
191{
192 return (n + k - 1) / k * k;
193}
194constexpr auto
195uround_up(usize n, usize k) noexcept -> usize
196{
197 return (n + k - 1) / k * k;
198}
199
200inline auto
201bytes_to_prev_aligned(void* ptr, usize align) noexcept -> isize
202{
203 using UPtr = std::uintptr_t;
204
205 UPtr mask = align - 1;
206 UPtr iptr = UPtr(ptr);
207 UPtr aligned_ptr = iptr & ~mask;
208 return isize(aligned_ptr - iptr);
209}
210inline auto
211bytes_to_next_aligned(void* ptr, usize align) noexcept -> isize
212{
213 using UPtr = std::uintptr_t;
214
215 UPtr mask = align - 1;
216 UPtr iptr = UPtr(ptr);
217 UPtr aligned_ptr = (iptr + mask) & ~mask;
218 return isize(aligned_ptr - iptr);
219}
220
221inline auto
222next_aligned(void* ptr, usize align) noexcept -> void*
223{
224 using BytePtr = unsigned char*;
225 using VoidPtr = void*;
226 return VoidPtr(BytePtr(ptr) + detail::bytes_to_next_aligned(ptr, align));
227}
228inline auto
229prev_aligned(void* ptr, usize align) noexcept -> void*
230{
231 using BytePtr = unsigned char*;
232 using VoidPtr = void*;
233 return VoidPtr(BytePtr(ptr) + detail::bytes_to_prev_aligned(ptr, align));
234}
235
236} // namespace detail
237
238enum struct Layout : unsigned char
239{
242};
243
246
247constexpr auto
249{
250 return Layout(1 - u32(l));
251}
252constexpr auto
254{
255 return l == colmajor ? Eigen::ColMajor : Eigen::RowMajor;
256}
257constexpr auto
259{
260 return (unsigned(l) & Eigen::RowMajorBit) == Eigen::RowMajor ? rowmajor
261 : colmajor;
262}
263
264static_assert(to_eigen_layout(from_eigen_layout(Eigen::ColMajor)) ==
265 Eigen::ColMajor,
266 ".");
267static_assert(to_eigen_layout(from_eigen_layout(Eigen::RowMajor)) ==
268 Eigen::RowMajor,
269 ".");
270
271namespace detail {
272template<Layout L>
274
275template<>
277{
278 template<typename T>
279 VEG_INLINE static constexpr auto offset(T* ptr,
280 isize row,
281 isize col,
282 isize outer_stride) noexcept -> T*
283 {
284 return ptr + (usize(row) + usize(col) * usize(outer_stride));
285 }
286
287 using NextRowStride = Eigen::Stride<0, 0>;
288 using NextColStride = Eigen::InnerStride<Eigen::Dynamic>;
289 VEG_INLINE static auto next_row_stride(isize outer_stride) noexcept
291 {
292 (void)outer_stride;
293 return NextRowStride{};
294 }
295 VEG_INLINE static auto next_col_stride(isize outer_stride) noexcept
297 {
298 return NextColStride /* NOLINT(modernize-return-braced-init-list) */ (
299 outer_stride);
300 }
301
302 template<typename T>
304 isize dim,
305 isize outer_stride)
306 {
307 (void)ptr, (void)dim, (void)outer_stride;
308 }
309};
310
311template<>
313{
314 template<typename T>
315 VEG_INLINE static constexpr auto offset(T* ptr,
316 isize row,
317 isize col,
318 isize outer_stride) noexcept -> T*
319 {
320 return ptr + (usize(col) + usize(row) * usize(outer_stride));
321 }
322
323 using NextColStride = Eigen::Stride<0, 0>;
324 using NextRowStride = Eigen::InnerStride<Eigen::Dynamic>;
325 VEG_INLINE static auto next_col_stride(isize outer_stride) noexcept
327 {
328 (void)outer_stride;
329 return NextColStride{};
330 }
331 VEG_INLINE static auto next_row_stride(isize outer_stride) noexcept
333 {
334 return NextRowStride /* NOLINT(modernize-return-braced-init-list) */ (
335 outer_stride);
336 }
337
338 template<typename T>
340 isize dim,
341 isize outer_stride)
342 {
343 Eigen::Map< //
344 Eigen::Matrix< //
345 T, //
346 Eigen::Dynamic, //
347 Eigen::Dynamic //
348 >, //
349 Eigen::Unaligned, //
350 Eigen::OuterStride<Eigen::Dynamic> //
351 >{
352 ptr,
353 dim,
354 dim,
355 Eigen::OuterStride<Eigen::Dynamic>(outer_stride),
356 }
357 .transposeInPlace();
358 }
359};
360} // namespace detail
361
362namespace detail {
363template<typename T>
364struct unlref
365{
366 using Type = T;
367};
368template<typename T>
369struct unlref<T&>
370{
371 using Type = T;
372};
373
374template<typename T>
375auto
376is_eigen_matrix_base_impl(Eigen::MatrixBase<T> const volatile*)
378auto
379is_eigen_matrix_base_impl(void const volatile*)
381
382template<typename T>
383auto
384is_eigen_owning_matrix_base_impl(Eigen::PlainObjectBase<T> const volatile*)
386auto
389
390template<typename... Ts>
391using Void = void;
392
393template<typename Mat, typename T>
394using DataExpr = decltype(static_cast<T*>(VEG_DECLVAL(Mat&).data()));
395
396template<typename Dummy,
397 typename Fallback,
398 template<typename...> class F,
399 typename... Ts>
404
405template<typename Fallback, template<typename...> class F, typename... Ts>
406struct DetectedImpl<Void<F<Ts...>>, Fallback, F, Ts...>
408{
409 using Type = F<Ts...>;
410};
411
412template<typename Fallback, template<typename...> class F, typename... Ts>
413using Detected = typename DetectedImpl<void, Fallback, F, Ts...>::Type;
414
415template<typename T>
417 proxsuite::linalg::veg::meta::constant<isize, isize(T::ColsAtCompileTime)>;
418template<typename T>
420 proxsuite::linalg::veg::meta::constant<isize, isize(T::RowsAtCompileTime)>;
421template<typename T>
424 isize(T::InnerStrideAtCompileTime)>;
425template<typename T>
426using LayoutImpl = proxsuite::linalg::veg::meta::
427 constant<Layout, (bool(T::IsRowMajor) ? rowmajor : colmajor)>;
428
429template<typename T, Layout L>
430using EigenMatMap = Eigen::Map< //
431 Eigen::Matrix< //
432 T, //
433 Eigen::Dynamic, //
434 Eigen::Dynamic, //
435 (L == colmajor) //
436 ? Eigen::ColMajor //
437 : Eigen::RowMajor //
438 > const, //
439 Eigen::Unaligned, //
440 Eigen::OuterStride<Eigen::Dynamic> //
441 >;
442template<typename T, Layout L>
443using EigenMatMapMut = Eigen::Map< //
444 Eigen::Matrix< //
445 T, //
446 Eigen::Dynamic, //
447 Eigen::Dynamic, //
448 (L == colmajor) //
449 ? Eigen::ColMajor //
450 : Eigen::RowMajor //
451 >, //
452 Eigen::Unaligned, //
453 Eigen::OuterStride<Eigen::Dynamic> //
454 >;
455
456template<typename T, typename Stride>
457using EigenVecMap = Eigen::Map< //
458 Eigen::Matrix< //
459 T, //
460 Eigen::Dynamic, //
461 1 //
462 > const, //
463 Eigen::Unaligned, //
464 Stride //
465 >;
466template<typename T, typename Stride>
467using EigenVecMapMut = Eigen::Map< //
468 Eigen::Matrix< //
469 T, //
470 Eigen::Dynamic, //
471 1 //
472 >, //
473 Eigen::Unaligned, //
474 Stride //
475 >;
476
477template<typename T, Layout L>
479template<typename T, Layout L>
481template<typename T, Layout L>
483template<typename T, Layout L>
485
486template<typename T>
488template<typename T>
490
491} // namespace detail
492
493template<typename T>
495
496namespace eigen {
497template<typename T>
501 T>;
502template<typename T>
506 T>;
507template<typename T>
511 T>;
512template<typename T>
514 detail::Detected<proxsuite::linalg::veg::meta::
515 constant<Layout, Layout(static_cast<unsigned char>(-1))>,
517 T>;
518} // namespace eigen
519
520namespace concepts {
521VEG_DEF_CONCEPT(typename T, rvalue_ref, std::is_rvalue_reference<T>::value);
522VEG_DEF_CONCEPT(typename T, lvalue_ref, std::is_lvalue_reference<T>::value);
523VEG_DEF_CONCEPT((template<typename...> class F, typename... Ts),
524 detected,
526
527namespace aux {
528VEG_DEF_CONCEPT((typename Mat, typename T),
529 has_data_expr,
530 LDLT_CONCEPT(detected<detail::DataExpr, Mat, T>));
531
532VEG_DEF_CONCEPT((typename Mat),
533 matrix_base,
535 static_cast<Mat*>(nullptr)))::value);
536
537VEG_DEF_CONCEPT((typename Mat),
538 is_plain_object_base,
540 static_cast<Mat*>(nullptr)))::value);
541
542VEG_DEF_CONCEPT((typename Mat),
543 tmp_matrix,
544 (LDLT_CONCEPT(aux::is_plain_object_base<unref<Mat>>) &&
545 !LDLT_CONCEPT(lvalue_ref<Mat>)));
546} // namespace aux
547
548VEG_DEF_CONCEPT((typename Mat, typename T),
549 eigen_view,
550 (LDLT_CONCEPT(aux::matrix_base<unref<Mat>>) &&
551 LDLT_CONCEPT(aux::has_data_expr<Mat, T const>)));
552
553VEG_DEF_CONCEPT((typename Mat, typename T),
554 eigen_view_mut,
555 (LDLT_CONCEPT(aux::matrix_base<unref<Mat>>) &&
556 LDLT_CONCEPT(aux::has_data_expr<Mat, T>) &&
557 !LDLT_CONCEPT(aux::tmp_matrix<Mat>)));
558
559VEG_DEF_CONCEPT((typename Mat, typename T),
560 eigen_strided_vector_view,
561 (LDLT_CONCEPT(eigen_view<Mat, T>) &&
562 (eigen::CompTimeCols<unref<Mat>>::value == 1)));
563
564VEG_DEF_CONCEPT((typename Mat, typename T),
565 eigen_strided_vector_view_mut,
566 (LDLT_CONCEPT(eigen_view_mut<Mat, T>) &&
567 (eigen::CompTimeCols<unref<Mat>>::value == 1)));
568
569VEG_DEF_CONCEPT((typename Mat, typename T),
570 eigen_vector_view,
571 (LDLT_CONCEPT(eigen_strided_vector_view<Mat, T>) &&
572 (eigen::CompTimeInnerStride<unref<Mat>>::value == 1)));
573
574VEG_DEF_CONCEPT((typename Mat, typename T),
575 eigen_vector_view_mut,
576 (LDLT_CONCEPT(eigen_strided_vector_view_mut<Mat, T>) &&
577 (eigen::CompTimeInnerStride<unref<Mat>>::value == 1)));
578} // namespace concepts
579
580inline namespace tags {
581VEG_TAG(from_ptr_size, FromPtrSize);
582VEG_TAG(from_ptr_size_stride, FromPtrSizeStride);
583VEG_TAG(from_ptr_rows_cols_stride, FromPtrRowsColsStride);
584VEG_TAG(from_eigen, FromEigen);
585} // namespace tags
586
587template<typename T>
589{
590 T const* data;
592
594 VectorView(FromPtrSize /*tag*/, T const* _data, isize _dim) noexcept
595 : data(_data)
596 , dim(_dim)
597 {
598 }
599
600 VEG_TEMPLATE(typename Vec,
601 requires(LDLT_CONCEPT(eigen_vector_view<Vec, T>)),
603 (/*tag*/, FromEigen),
604 (vec, Vec const&)) noexcept
605 : data(vec.data())
606 , dim(vec.rows())
607 {
608 }
609
610 VEG_INLINE auto ptr(isize index) const noexcept -> T const*
611 {
612 return data + index;
613 }
614 VEG_INLINE auto operator()(isize index) const noexcept -> T const&
615 {
616 return *ptr(index);
617 }
618 VEG_INLINE auto segment(isize i, isize size) const noexcept -> VectorView
619 {
620 return {
621 from_ptr_size,
622 data + i,
623 size,
624 };
625 }
626 VEG_INLINE auto to_eigen() const -> detail::VecMap<T>
627 {
628 return detail::VecMap<T>(data, Eigen::Index(dim));
629 }
630};
631
632template<typename T>
634{
637
639 VectorViewMut(FromPtrSize /*tag*/, T* _data, isize _dim) noexcept
640 : data(_data)
641 , dim(_dim)
642 {
643 }
644
645 VEG_TEMPLATE(typename Vec,
646 requires(LDLT_CONCEPT(eigen_vector_view_mut<Vec, T>)),
648 (/*tag*/, FromEigen),
649 (vec, Vec&&)) noexcept
650 : data(vec.data())
651 , dim(vec.rows())
652 {
653 }
654
655 VEG_INLINE auto as_const() const noexcept -> VectorView<T>
656 {
657 return {
658 from_ptr_size,
659 data,
660 dim,
661 };
662 }
663 VEG_INLINE auto ptr(isize index) const noexcept -> T* { return data + index; }
664 VEG_INLINE auto operator()(isize index) const noexcept -> T&
665 {
666 return *ptr(index);
667 }
668 VEG_INLINE auto segment(isize i, isize size) const noexcept -> VectorViewMut
669 {
670 return {
671 from_ptr_size,
672 data + i,
673 size,
674 };
675 }
676 VEG_INLINE auto to_eigen() const -> detail::VecMapMut<T>
677 {
678 return detail::VecMapMut<T>(data, Eigen::Index(dim));
679 }
680};
681
682template<typename T>
684{
685 T const* data;
688
690 StridedVectorView(FromPtrSizeStride /*tag*/,
691 T const* _data,
692 isize _dim,
693 isize _stride) noexcept
694 : data(_data)
695 , dim(_dim)
696 , stride(_stride)
697 {
698 }
699
700 VEG_TEMPLATE(typename Vec,
701 requires(LDLT_CONCEPT(eigen_strided_vector_view<Vec, T>)),
703 (/*tag*/, FromEigen),
704 (vec, Vec const&)) noexcept
705 : data(vec.data())
706 , dim(vec.rows())
707 , stride(vec.innerStride())
708 {
709 }
710
711 VEG_INLINE auto ptr(isize index) const noexcept -> T const*
712 {
713 return data + stride * index;
714 }
715 VEG_INLINE auto operator()(isize index) const noexcept -> T const&
716 {
717 return *ptr(index);
718 }
719 VEG_INLINE auto segment(isize i, isize size) const noexcept
721 {
722 return {
723 from_ptr_size_stride,
724 data + stride * i,
725 size,
726 stride,
727 };
728 }
729 VEG_INLINE auto to_eigen() const
730 -> detail::EigenVecMap<T, Eigen::InnerStride<Eigen::Dynamic>>
731 {
733 data,
734 Eigen::Index(dim),
735 Eigen::Index(1),
736 Eigen::InnerStride<Eigen::Dynamic>(Eigen::Index(stride)));
737 }
738};
739
740template<typename T>
742{
746
748 StridedVectorViewMut(FromPtrSizeStride /*tag*/,
749 T* _data,
750 isize _dim,
751 isize _stride) noexcept
752 : data(_data)
753 , dim(_dim)
754 , stride(_stride)
755 {
756 }
757
758 VEG_TEMPLATE(typename Vec,
759 requires(LDLT_CONCEPT(eigen_strided_vector_view_mut<Vec, T>)),
761 (/*tag*/, FromEigen),
762 (vec, Vec&&)) noexcept
763 : data(vec.data())
764 , dim(vec.rows())
765 , stride(vec.innerStride())
766 {
767 }
768
769 VEG_INLINE auto as_const() const noexcept -> StridedVectorView<T>
770 {
771 return {
772 from_ptr_size_stride,
773 data,
774 dim,
775 stride,
776 };
777 }
778 VEG_INLINE auto ptr(isize index) const noexcept -> T*
779 {
780 return data + stride * index;
781 }
782 VEG_INLINE auto operator()(isize index) const noexcept -> T&
783 {
784 return *ptr(index);
785 }
786 VEG_INLINE auto segment(isize i, isize size) const noexcept
788 {
789 return {
790 from_ptr_size_stride,
791 data + stride * i,
792 size,
793 stride,
794 };
795 }
796 VEG_INLINE auto to_eigen() const
797 -> detail::EigenVecMapMut<T, Eigen::InnerStride<Eigen::Dynamic>>
798 {
800 data,
801 Eigen::Index(dim),
802 Eigen::Index(1),
803 Eigen::InnerStride<Eigen::Dynamic>(Eigen::Index(stride)));
804 }
805};
806
807template<typename T, Layout L>
809{
810 T const* data;
814
815 VEG_INLINE MatrixView(FromPtrRowsColsStride /*tag*/,
816 T const* _data,
817 isize _rows,
818 isize _cols,
819 isize _outer_stride) noexcept
820 : data(_data)
821 , rows(_rows)
822 , cols(_cols)
823 , outer_stride(_outer_stride)
824 {
825 }
826
827 VEG_TEMPLATE(typename Mat,
828 requires(LDLT_CONCEPT(eigen_view<Mat, T>) &&
829 eigen::GetLayout<unref<Mat>>::value == L),
831 (/*tag*/, FromEigen),
832 (mat, Mat const&)) noexcept
833 : data(mat.data())
834 , rows(mat.rows())
835 , cols(mat.cols())
836 , outer_stride(mat.outerStride())
837 {
838 }
839
840 VEG_INLINE auto ptr(isize row, isize col) const noexcept -> T const*
841 {
843 }
844 VEG_INLINE auto operator()(isize row, isize col) const noexcept -> T const&
845 {
846 return *ptr(row, col);
847 }
849 isize col,
850 isize nrows,
851 isize ncols) const noexcept -> MatrixView
852 {
853 return {
854 from_ptr_rows_cols_stride,
856 nrows,
857 ncols,
859 };
860 }
861
862private:
863 VEG_INLINE auto col_impl(
865 isize c) const noexcept -> VectorView<T>
866 {
867 return {
868 from_ptr_size,
869 data + c * outer_stride,
870 rows,
871 };
872 }
873 VEG_INLINE auto col_impl(
875 isize c) const noexcept -> StridedVectorView<T>
876 {
877 return {
878 from_ptr_size_stride,
879 data + c,
880 rows,
882 };
883 }
884
885public:
886 VEG_INLINE auto col(isize c) const noexcept -> proxsuite::linalg::veg::meta::
888 {
890 }
891 VEG_INLINE auto row(isize r) const noexcept -> proxsuite::linalg::veg::meta::
893 {
894 return trans().col(r);
895 }
896 VEG_INLINE auto trans() const noexcept
897 -> MatrixView<T, proxqp::flip_layout(L)>
898 {
899 return {
900 from_ptr_rows_cols_stride, data, cols, rows, outer_stride,
901 };
902 }
903 VEG_INLINE auto to_eigen() const noexcept -> detail::EigenMatMap<T, L>
904 {
906 data,
907 Eigen::Index(rows),
908 Eigen::Index(cols),
909 Eigen::OuterStride<Eigen::Dynamic>(Eigen::Index(outer_stride)));
910 }
911};
912
913template<typename T, Layout L>
915{
920
921 VEG_INLINE MatrixViewMut(FromPtrRowsColsStride /*tag*/,
922 T* _data,
923 isize _rows,
924 isize _cols,
925 isize _outer_stride) noexcept
926 : data(_data)
927 , rows(_rows)
928 , cols(_cols)
929 , outer_stride(_outer_stride)
930 {
931 }
932
933 VEG_TEMPLATE(typename Mat,
934 requires(LDLT_CONCEPT(eigen_view<Mat, T>) &&
935 eigen::GetLayout<unref<Mat>>::value == L),
937 (/*tag*/, FromEigen),
938 (mat, Mat&&)) noexcept
939 : data(mat.data())
940 , rows(mat.rows())
941 , cols(mat.cols())
942 , outer_stride(mat.outerStride())
943 {
944 }
945
946 VEG_INLINE auto ptr(isize row, isize col) const noexcept -> T*
947 {
949 }
950 VEG_INLINE auto operator()(isize row, isize col) const noexcept -> T&
951 {
952 return *ptr(row, col);
953 }
955 isize col,
956 isize nrows,
957 isize ncols) const noexcept -> MatrixViewMut
958 {
959 return {
960 from_ptr_rows_cols_stride,
962 nrows,
963 ncols,
965 };
966 }
967
968private:
969 VEG_INLINE auto col_impl(
971 isize c) const noexcept -> VectorViewMut<T>
972 {
973 return {
974 from_ptr_size,
975 data + c * outer_stride,
976 rows,
977 };
978 }
979 VEG_INLINE auto col_impl(
981 isize c) const noexcept -> StridedVectorViewMut<T>
982 {
983 return {
984 from_ptr_size_stride,
985 data + c,
986 rows,
988 };
989 }
990
991public:
992 VEG_INLINE auto col(isize c) const noexcept -> proxsuite::linalg::veg::meta::
994 {
996 }
997 VEG_INLINE auto row(isize r) const noexcept -> proxsuite::linalg::veg::meta::
999 {
1000 return trans().col(r);
1001 }
1002 VEG_INLINE auto trans() const noexcept
1004 {
1005 return {
1006 from_ptr_rows_cols_stride, data, cols, rows, outer_stride,
1007 };
1008 }
1009 VEG_INLINE auto to_eigen() const noexcept -> detail::EigenMatMapMut<T, L>
1010 {
1012 data,
1013 Eigen::Index(rows),
1014 Eigen::Index(cols),
1015 Eigen::OuterStride<Eigen::Dynamic>(Eigen::Index(outer_stride)));
1016 }
1017 VEG_INLINE auto as_const() const noexcept -> MatrixView<T, L>
1018 {
1019 return {
1020 from_ptr_rows_cols_stride, data, rows, cols, outer_stride,
1021 };
1022 }
1023};
1024
1025template<typename T>
1027{
1028private:
1030
1031public:
1032 explicit LdltView(MatrixView<T, colmajor> ld) noexcept
1033 : ld(ld)
1034 {
1035 VEG_DEBUG_ASSERT(ld.rows == ld.cols);
1036 }
1037
1038 VEG_INLINE auto l() const noexcept -> MatrixView<T, colmajor> { return ld; }
1039 VEG_INLINE auto d() const noexcept -> StridedVectorView<T>
1040 {
1041 return { from_ptr_size_stride, ld.data, ld.rows, ld.outer_stride + 1 };
1042 }
1043
1045 {
1046 return LdltView{ ld.block(0, 0, k, k) };
1047 }
1049 {
1050 isize n = ld.rows;
1051 return LdltView{ ld.block(n - k, n - k, k, k) };
1052 }
1053};
1054template<typename T>
1056{
1057private:
1059
1060public:
1062 : ld(ld)
1063 {
1064 VEG_DEBUG_ASSERT(ld.rows == ld.cols);
1065 }
1066
1067 VEG_INLINE auto l() const noexcept -> MatrixView<T, colmajor>
1068 {
1069 return ld.as_const();
1070 }
1071 VEG_INLINE auto l_mut() const noexcept -> MatrixViewMut<T, colmajor>
1072 {
1073 return ld;
1074 }
1075 VEG_INLINE auto d() const noexcept -> StridedVectorView<T>
1076 {
1077 return { from_ptr_size_stride, ld.data, ld.rows, ld.outer_stride + 1 };
1078 }
1079 VEG_INLINE auto d_mut() const noexcept -> StridedVectorViewMut<T>
1080 {
1081 return { from_ptr_size_stride, ld.data, ld.rows, ld.outer_stride + 1 };
1082 }
1083
1084 VEG_INLINE auto as_const() const noexcept -> LdltView<T>
1085 {
1086 return LdltView<T>{ ld.as_const() };
1087 }
1088
1090 {
1091 return LdltViewMut{ ld.block(0, 0, k, k) };
1092 }
1094 {
1095 isize n = ld.rows;
1096 return LdltViewMut{ ld.block(n - k, n - k, k, k) };
1097 }
1098};
1099
1100namespace detail {
1101template<typename T>
1102void
1106 T factor)
1107{
1108
1109 if ((dst.cols == 0) || (dst.rows == 0) || (lhs.cols == 0)) {
1110 return;
1111 }
1112
1113#if !EIGEN_VERSION_AT_LEAST(3, 3, 8)
1114#define LAZY_PRODUCT(a, b) a.lazyProduct(b)
1115#else
1116#define LAZY_PRODUCT(a, b) a.operator*(b)
1117#endif
1118
1119 if (dst.cols == 1 && dst.rows == 1) {
1120 // dot
1121 auto rhs_col = rhs.col(0);
1122 auto lhs_row = lhs.row(0);
1123 auto lhs_as_col = lhs.col(0);
1124 lhs_as_col.dim = lhs_row.dim;
1125 if (lhs_row.stride == 1) {
1126 dst(0, 0) += factor * lhs_as_col.to_eigen().dot(rhs_col.to_eigen());
1127 } else {
1128 dst(0, 0) += factor * lhs_row.to_eigen().dot(rhs_col.to_eigen());
1129 }
1130 } else if (dst.cols == 1) {
1131 // gemv
1132 auto rhs_col = rhs.col(0);
1133 auto dst_col = dst.col(0);
1134 dst_col.to_eigen().noalias().operator+=(
1135 factor * LAZY_PRODUCT(lhs.to_eigen(), rhs_col.to_eigen()));
1136 }
1137
1138#if !EIGEN_VERSION_AT_LEAST(3, 3, 8)
1139 else if ((dst.rows < 20) && (dst.cols < 20) && (rhs.rows < 20)) {
1140 // gemm
1141 // workaround for eigen 3.3.7 bug:
1142 // https://gitlab.com/libeigen/eigen/-/issues/1562
1143 using Stride = Eigen::OuterStride<Eigen::Dynamic>;
1144 using Mat =
1145 Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::ColMajor, 20, 20>;
1146 using MapMut = Eigen::Map<Mat, Eigen::Unaligned, Stride>;
1147 using Map = Eigen::Map<Mat const, Eigen::Unaligned, Stride>;
1148
1149 MapMut(dst.data, dst.rows, dst.cols, Stride(dst.outer_stride))
1150 .noalias()
1151 .operator+=(
1152 factor *
1154 Map(lhs.data, lhs.rows, lhs.cols, Stride(lhs.outer_stride)),
1155 Map(rhs.data, rhs.rows, rhs.cols, Stride(rhs.outer_stride))));
1156 }
1157#endif
1158
1159 else {
1160 // gemm
1161 dst.to_eigen().noalias().operator+=(
1162 factor * LAZY_PRODUCT(lhs.to_eigen(), rhs.to_eigen()));
1163 }
1164
1165#undef LAZY_PRODUCT
1166}
1167
1168template<typename T>
1169void
1172 VectorView<T> rhs,
1173 T factor)
1174{
1176 {
1177 from_ptr_rows_cols_stride,
1178 dst.data,
1179 dst.dim,
1180 1,
1181 0,
1182 },
1183 lhs,
1184 {
1185 from_ptr_rows_cols_stride,
1186 rhs.data,
1187 rhs.dim,
1188 1,
1189 0,
1190 },
1191 VEG_FWD(factor));
1192}
1193
1194template<typename T>
1195auto
1197{
1198 auto out = T(0);
1200 {
1201 from_ptr_rows_cols_stride,
1202 std::addressof(out),
1203 1,
1204 1,
1205 0,
1206 },
1207 {
1208 from_ptr_rows_cols_stride,
1209 lhs.data,
1210 1,
1211 lhs.dim,
1212 lhs.stride,
1213 },
1214 {
1215 from_ptr_rows_cols_stride,
1216 rhs.data,
1217 rhs.dim,
1218 1,
1219 0,
1220 },
1221 1);
1222 return out;
1223}
1224template<typename T>
1225void
1229{
1230 out.to_eigen() = lhs.to_eigen().cwiseProduct(rhs.to_eigen());
1231}
1232template<typename T>
1233void
1235{
1236 out.to_eigen() = in.to_eigen().operator*(factor);
1237}
1238
1239template<typename T>
1240void
1243{
1244 if (rhs.cols == 1) {
1245 tr.to_eigen()
1246 .transpose()
1247 .template triangularView<Eigen::UnitUpper>()
1248 .template solveInPlace<Eigen::OnTheRight>(rhs.col(0).to_eigen());
1249 } else {
1250 tr.to_eigen()
1251 .transpose()
1252 .template triangularView<Eigen::UnitUpper>()
1253 .template solveInPlace<Eigen::OnTheRight>(rhs.to_eigen());
1254 }
1255}
1256
1257template<typename T>
1258void
1262{
1263 if (out.cols == 1) {
1264 out.col(0).to_eigen() =
1265 in.col(0).to_eigen().operator*(d.to_eigen().asDiagonal().inverse());
1266 } else {
1267 out.to_eigen() =
1268 in.to_eigen().operator*(d.to_eigen().asDiagonal().inverse());
1269 }
1270}
1271template<typename T>
1272void
1276{
1277 if (out.cols == 1) {
1278 out.col(0).to_eigen() =
1279 in.col(0).to_eigen().operator*(d.to_eigen().asDiagonal());
1280 } else {
1281 out.to_eigen() = in.to_eigen().operator*(d.to_eigen().asDiagonal());
1282 }
1283}
1284
1285template<typename T>
1286void
1290{
1291 if (lhs.cols == 1) {
1292 out.to_eigen().template triangularView<Eigen::Lower>().operator-=(
1293 lhs.col(0).to_eigen().operator*(
1294 Eigen::Map<Eigen::Matrix<T, 1, Eigen::Dynamic> const>(
1295 rhs.data, 1, rhs.cols)));
1296 } else {
1297 out.to_eigen().template triangularView<Eigen::Lower>().operator-=(
1298 lhs.to_eigen().operator*(rhs.to_eigen()));
1299 }
1300}
1301
1302} // namespace detail
1303} // namespace proxqp
1304} // namespace proxsuite
1305
1306namespace proxsuite {
1307namespace proxqp {
1308
1309namespace dense {
1310
1312{
1317 auto operator=(EigenAllowAlloc const&) -> EigenAllowAlloc& = delete;
1318
1319#if defined(EIGEN_RUNTIME_NO_MALLOC)
1320 EigenAllowAlloc() noexcept
1321 : alloc_was_allowed(Eigen::internal::is_malloc_allowed())
1322 {
1323 Eigen::internal::set_is_malloc_allowed(true);
1324 }
1325 ~EigenAllowAlloc() noexcept
1326 {
1327 Eigen::internal::set_is_malloc_allowed(alloc_was_allowed);
1328 }
1329#else
1330 EigenAllowAlloc() = default;
1331#endif
1332};
1333
1334template<typename T>
1348
1349template<typename Scalar>
1366
1367template<typename T>
1369{
1370 static constexpr Layout layout = rowmajor;
1371
1374
1379
1380 VEG_INLINE constexpr auto as_const() const noexcept -> QpView<T>
1381 {
1382 return {
1383 H.as_const(), g.as_const(), A.as_const(),
1384 b.as_const(), C.as_const(), d.as_const(),
1385 };
1386 }
1387};
1388
1389template<typename Scalar>
1391{
1392 static constexpr Layout layout = rowmajor;
1393
1396
1405
1406 VEG_INLINE constexpr auto as_const() const noexcept -> QpViewBox<Scalar>
1407 {
1408 return { H.as_const(), g.as_const(), A.as_const(), b.as_const(),
1409 C.as_const(), u.as_const(), l.as_const(), I.as_const(),
1410 u_box.as_const(), l_box.as_const() };
1411 }
1412};
1413
1414namespace nb {
1415struct pow
1416{
1417 template<typename T>
1418 auto operator()(T x, T y) const -> T
1419 {
1420 using std::pow;
1421 return pow(x, y);
1422 }
1423};
1425{
1426 template<typename D>
1427 auto operator()(Eigen::MatrixBase<D> const& mat) const -> typename D::Scalar
1428 {
1429 if (mat.rows() == 0 || mat.cols() == 0) {
1430 return typename D::Scalar(0);
1431 } else {
1432 return mat.template lpNorm<Eigen::Infinity>();
1433 }
1434 }
1435};
1436struct sqrt
1437{
1438 template<typename T>
1439 auto operator()(T x) const -> T
1440 {
1441 using std::sqrt;
1442 return sqrt(x);
1443 }
1444};
1445struct fabs
1446{
1447 template<typename T>
1448 auto operator()(T x) const -> T
1449 {
1450 using std::fabs;
1451 return fabs(x);
1452 }
1453};
1454} // namespace nb
1458VEG_NIEBLOID(infty_norm);
1459} // namespace dense
1460} // namespace proxqp
1461} // namespace proxsuite
1462
1463#endif /* end of include guard PROXSUITE_PROXQP_DENSE_VIEWS_HPP */
#define VEG_DEBUG_ASSERT(...)
Definition assert.hpp:38
#define LAZY_PRODUCT(a, b)
#define LDLT_CONCEPT(...)
Definition views.hpp:16
#define VEG_TAG(Name, Type)
Definition macros.hpp:550
#define VEG_NIEBLOID(Name)
Definition macros.hpp:545
#define VEG_INLINE
Definition macros.hpp:118
#define VEG_FWD(X)
Definition macros.hpp:569
#define VEG_DECLVAL(...)
Definition macros.hpp:131
#define VEG_DEF_CONCEPT(Tpl, Name,...)
Definition macros.hpp:321
bool_constant< false > false_type
Definition macros.hpp:903
bool_constant< true > true_type
Definition macros.hpp:902
decltype(sizeof(0)) usize
Definition macros.hpp:702
typename _detail::pack_ith_elem< I >::template Type< Ts... > ith
std::uint64_t u64
Definition typedefs.hpp:46
std::uint32_t u32
Definition typedefs.hpp:48
_detail::_meta::make_signed< usize >::Type isize
Definition typedefs.hpp:43
proxsuite::linalg::veg::meta::constant< isize, isize(T::RowsAtCompileTime)> CompTimeRowsImpl
Definition views.hpp:419
auto prev_aligned(void *ptr, usize align) noexcept -> void *
Definition views.hpp:229
void apply_diag_on_right(MatrixViewMut< T, colmajor > out, StridedVectorView< T > d, MatrixView< T, colmajor > in)
Definition views.hpp:1273
EigenVecMapMut< T, typename ElementAccess< L >::NextRowStride > ColToVecMut
Definition views.hpp:482
constexpr auto cx_min_list(T const (&arr)[N]) noexcept -> T
Definition views.hpp:129
auto is_eigen_matrix_base_impl(Eigen::MatrixBase< T > const volatile *) -> proxsuite::linalg::veg::meta::true_type
auto is_eigen_owning_matrix_base_impl(Eigen::PlainObjectBase< T > const volatile *) -> proxsuite::linalg::veg::meta::true_type
Eigen::Map< Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic,(L==colmajor) ? Eigen::ColMajor :Eigen::RowMajor > const, Eigen::Unaligned, Eigen::OuterStride< Eigen::Dynamic > > EigenMatMap
Definition views.hpp:430
EigenVecMap< T, Eigen::Stride< 0, 0 > > VecMap
Definition views.hpp:487
auto bytes_to_prev_aligned(void *ptr, usize align) noexcept -> isize
Definition views.hpp:201
auto next_aligned(void *ptr, usize align) noexcept -> void *
Definition views.hpp:222
EigenVecMapMut< T, Eigen::Stride< 0, 0 > > VecMapMut
Definition views.hpp:489
void noalias_mul_add(MatrixViewMut< T, colmajor > dst, MatrixView< T, colmajor > lhs, MatrixView< T, colmajor > rhs, T factor)
Definition views.hpp:1103
auto bytes_to_next_aligned(void *ptr, usize align) noexcept -> isize
Definition views.hpp:211
constexpr auto min_list_impl(T init, T const *arr, usize n) noexcept -> T
Definition views.hpp:121
EigenVecMap< T, typename ElementAccess< L >::NextRowStride > ColToVec
Definition views.hpp:478
auto dot(StridedVectorView< T > lhs, VectorView< T > rhs) -> T
Definition views.hpp:1196
proxsuite::linalg::veg::meta::constant< isize, isize(T::ColsAtCompileTime)> CompTimeColsImpl
Definition views.hpp:416
decltype(static_cast< T * >(VEG_DECLVAL(Mat &).data())) DataExpr
Definition views.hpp:394
Eigen::Map< Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic,(L==colmajor) ? Eigen::ColMajor :Eigen::RowMajor >, Eigen::Unaligned, Eigen::OuterStride< Eigen::Dynamic > > EigenMatMapMut
Definition views.hpp:443
void set_zero(T *dest, usize n)
Definition views.hpp:184
void apply_diag_inv_on_right(MatrixViewMut< T, colmajor > out, StridedVectorView< T > d, MatrixView< T, colmajor > in)
Definition views.hpp:1259
proxsuite::linalg::veg::meta::constant< isize, isize(T::InnerStrideAtCompileTime)> CompTimeInnerStrideImpl
Definition views.hpp:422
EigenVecMapMut< T, typename ElementAccess< L >::NextColStride > RowToVecMut
Definition views.hpp:484
EigenVecMap< T, typename ElementAccess< L >::NextColStride > RowToVec
Definition views.hpp:480
Eigen::Map< Eigen::Matrix< T, Eigen::Dynamic, 1 >, Eigen::Unaligned, Stride > EigenVecMapMut
Definition views.hpp:467
void assign_scalar_prod(VectorViewMut< T > out, T factor, VectorView< T > in)
Definition views.hpp:1234
proxsuite::linalg::veg::meta:: constant< Layout,(bool(T::IsRowMajor) ? rowmajor :colmajor)> LayoutImpl
Definition views.hpp:426
void trans_tr_unit_up_solve_in_place_on_right(MatrixView< T, colmajor > tr, MatrixViewMut< T, colmajor > rhs)
Definition views.hpp:1241
constexpr auto round_up(isize n, isize k) noexcept -> isize
Definition views.hpp:190
void assign_cwise_prod(VectorViewMut< T > out, StridedVectorView< T > lhs, StridedVectorView< T > rhs)
Definition views.hpp:1226
Eigen::Map< Eigen::Matrix< T, Eigen::Dynamic, 1 > const, Eigen::Unaligned, Stride > EigenVecMap
Definition views.hpp:457
void noalias_mul_add_vec(VectorViewMut< T > dst, MatrixView< T, colmajor > lhs, VectorView< T > rhs, T factor)
Definition views.hpp:1170
void noalias_mul_sub_tr_lo(MatrixViewMut< T, colmajor > out, MatrixView< T, colmajor > lhs, MatrixView< T, rowmajor > rhs)
Definition views.hpp:1287
typename DetectedImpl< void, Fallback, F, Ts... >::Type Detected
Definition views.hpp:413
constexpr auto uround_up(usize n, usize k) noexcept -> usize
Definition views.hpp:195
detail::Detected< proxsuite::linalg::veg::meta::constant< isize, 0 >, detail::CompTimeColsImpl, T > CompTimeCols
Definition views.hpp:498
detail::Detected< proxsuite::linalg::veg::meta::constant< isize, 0 >, detail::CompTimeInnerStrideImpl, T > CompTimeInnerStride
Definition views.hpp:508
detail::Detected< proxsuite::linalg::veg::meta::constant< isize, 0 >, detail::CompTimeRowsImpl, T > CompTimeRows
Definition views.hpp:503
detail::Detected< proxsuite::linalg::veg::meta:: constant< Layout, Layout(static_cast< unsigned char >(-1))>, detail::LayoutImpl, T > GetLayout
Definition views.hpp:513
constexpr auto flip_layout(Layout l) noexcept -> Layout
Definition views.hpp:248
constexpr auto from_eigen_layout(int l) -> Layout
Definition views.hpp:258
constexpr Layout colmajor
Definition views.hpp:244
constexpr auto to_eigen_layout(Layout l) -> int
Definition views.hpp:253
typename detail::unlref< T & >::Type unref
Definition views.hpp:494
constexpr Layout rowmajor
Definition views.hpp:245
VEG_INLINE auto l_mut() const noexcept -> MatrixViewMut< T, colmajor >
Definition views.hpp:1071
LdltViewMut(MatrixViewMut< T, colmajor > ld) noexcept
Definition views.hpp:1061
VEG_INLINE auto head(isize k) const -> LdltViewMut
Definition views.hpp:1089
VEG_INLINE auto d() const noexcept -> StridedVectorView< T >
Definition views.hpp:1075
VEG_INLINE auto d_mut() const noexcept -> StridedVectorViewMut< T >
Definition views.hpp:1079
VEG_INLINE auto as_const() const noexcept -> LdltView< T >
Definition views.hpp:1084
VEG_INLINE auto l() const noexcept -> MatrixView< T, colmajor >
Definition views.hpp:1067
VEG_INLINE auto tail(isize k) const -> LdltViewMut
Definition views.hpp:1093
VEG_INLINE auto d() const noexcept -> StridedVectorView< T >
Definition views.hpp:1039
VEG_INLINE auto head(isize k) const -> LdltView
Definition views.hpp:1044
LdltView(MatrixView< T, colmajor > ld) noexcept
Definition views.hpp:1032
VEG_INLINE auto l() const noexcept -> MatrixView< T, colmajor >
Definition views.hpp:1038
VEG_INLINE auto tail(isize k) const -> LdltView
Definition views.hpp:1048
VEG_INLINE auto row(isize r) const noexcept -> proxsuite::linalg::veg::meta::if_t<(L==rowmajor), VectorViewMut< T >, StridedVectorViewMut< T > >
Definition views.hpp:997
VEG_INLINE auto operator()(isize row, isize col) const noexcept -> T &
Definition views.hpp:950
VEG_TEMPLATE(typename Mat, requires(LDLT_CONCEPT(eigen_view< Mat, T >) &&eigen::GetLayout< unref< Mat > >::value==L), VEG_INLINE MatrixViewMut,(, FromEigen),(mat, Mat &&)) noexcept
Definition views.hpp:933
VEG_INLINE auto trans() const noexcept -> MatrixViewMut< T, proxqp::flip_layout(L)>
Definition views.hpp:1002
VEG_INLINE auto to_eigen() const noexcept -> detail::EigenMatMapMut< T, L >
Definition views.hpp:1009
VEG_INLINE auto as_const() const noexcept -> MatrixView< T, L >
Definition views.hpp:1017
VEG_INLINE auto block(isize row, isize col, isize nrows, isize ncols) const noexcept -> MatrixViewMut
Definition views.hpp:954
VEG_INLINE auto col(isize c) const noexcept -> proxsuite::linalg::veg::meta::if_t<(L==colmajor), VectorViewMut< T >, StridedVectorViewMut< T > >
Definition views.hpp:992
VEG_INLINE MatrixViewMut(FromPtrRowsColsStride, T *_data, isize _rows, isize _cols, isize _outer_stride) noexcept
Definition views.hpp:921
VEG_INLINE auto ptr(isize row, isize col) const noexcept -> T *
Definition views.hpp:946
VEG_INLINE auto block(isize row, isize col, isize nrows, isize ncols) const noexcept -> MatrixView
Definition views.hpp:848
VEG_INLINE auto operator()(isize row, isize col) const noexcept -> T const &
Definition views.hpp:844
VEG_INLINE auto to_eigen() const noexcept -> detail::EigenMatMap< T, L >
Definition views.hpp:903
VEG_INLINE MatrixView(FromPtrRowsColsStride, T const *_data, isize _rows, isize _cols, isize _outer_stride) noexcept
Definition views.hpp:815
VEG_INLINE auto col(isize c) const noexcept -> proxsuite::linalg::veg::meta::if_t<(L==colmajor), VectorView< T >, StridedVectorView< T > >
Definition views.hpp:886
VEG_INLINE auto row(isize r) const noexcept -> proxsuite::linalg::veg::meta::if_t<(L==rowmajor), VectorView< T >, StridedVectorView< T > >
Definition views.hpp:891
VEG_INLINE auto trans() const noexcept -> MatrixView< T, proxqp::flip_layout(L)>
Definition views.hpp:896
VEG_TEMPLATE(typename Mat, requires(LDLT_CONCEPT(eigen_view< Mat, T >) &&eigen::GetLayout< unref< Mat > >::value==L), VEG_INLINE MatrixView,(, FromEigen),(mat, Mat const &)) noexcept
Definition views.hpp:827
VEG_INLINE auto ptr(isize row, isize col) const noexcept -> T const *
Definition views.hpp:840
VEG_TEMPLATE(typename Vec, requires(LDLT_CONCEPT(eigen_strided_vector_view_mut< Vec, T >)), VEG_INLINE StridedVectorViewMut,(, FromEigen),(vec, Vec &&)) noexcept
Definition views.hpp:758
VEG_INLINE auto as_const() const noexcept -> StridedVectorView< T >
Definition views.hpp:769
VEG_INLINE auto operator()(isize index) const noexcept -> T &
Definition views.hpp:782
VEG_INLINE auto ptr(isize index) const noexcept -> T *
Definition views.hpp:778
VEG_INLINE StridedVectorViewMut(FromPtrSizeStride, T *_data, isize _dim, isize _stride) noexcept
Definition views.hpp:748
VEG_INLINE auto to_eigen() const -> detail::EigenVecMapMut< T, Eigen::InnerStride< Eigen::Dynamic > >
Definition views.hpp:796
VEG_INLINE auto segment(isize i, isize size) const noexcept -> StridedVectorViewMut
Definition views.hpp:786
VEG_INLINE auto operator()(isize index) const noexcept -> T const &
Definition views.hpp:715
VEG_INLINE auto segment(isize i, isize size) const noexcept -> StridedVectorView
Definition views.hpp:719
VEG_INLINE auto ptr(isize index) const noexcept -> T const *
Definition views.hpp:711
VEG_INLINE auto to_eigen() const -> detail::EigenVecMap< T, Eigen::InnerStride< Eigen::Dynamic > >
Definition views.hpp:729
VEG_TEMPLATE(typename Vec, requires(LDLT_CONCEPT(eigen_strided_vector_view< Vec, T >)), VEG_INLINE StridedVectorView,(, FromEigen),(vec, Vec const &)) noexcept
Definition views.hpp:700
VEG_INLINE StridedVectorView(FromPtrSizeStride, T const *_data, isize _dim, isize _stride) noexcept
Definition views.hpp:690
VEG_INLINE auto segment(isize i, isize size) const noexcept -> VectorViewMut
Definition views.hpp:668
VEG_TEMPLATE(typename Vec, requires(LDLT_CONCEPT(eigen_vector_view_mut< Vec, T >)), VEG_INLINE VectorViewMut,(, FromEigen),(vec, Vec &&)) noexcept
Definition views.hpp:645
VEG_INLINE auto operator()(isize index) const noexcept -> T &
Definition views.hpp:664
VEG_INLINE auto to_eigen() const -> detail::VecMapMut< T >
Definition views.hpp:676
VEG_INLINE auto ptr(isize index) const noexcept -> T *
Definition views.hpp:663
VEG_INLINE auto as_const() const noexcept -> VectorView< T >
Definition views.hpp:655
VEG_INLINE VectorViewMut(FromPtrSize, T *_data, isize _dim) noexcept
Definition views.hpp:639
VEG_INLINE auto ptr(isize index) const noexcept -> T const *
Definition views.hpp:610
VEG_INLINE VectorView(FromPtrSize, T const *_data, isize _dim) noexcept
Definition views.hpp:594
VEG_INLINE auto segment(isize i, isize size) const noexcept -> VectorView
Definition views.hpp:618
VEG_INLINE auto to_eigen() const -> detail::VecMap< T >
Definition views.hpp:626
VEG_INLINE auto operator()(isize index) const noexcept -> T const &
Definition views.hpp:614
VEG_TEMPLATE(typename Vec, requires(LDLT_CONCEPT(eigen_vector_view< Vec, T >)), VEG_INLINE VectorView,(, FromEigen),(vec, Vec const &)) noexcept
Definition views.hpp:600
EigenAllowAlloc(EigenAllowAlloc &&)=delete
auto operator=(EigenAllowAlloc const &) -> EigenAllowAlloc &=delete
EigenAllowAlloc(EigenAllowAlloc const &)=delete
auto operator=(EigenAllowAlloc &&) -> EigenAllowAlloc &=delete
MatrixViewMut< Scalar, layout > C
Definition views.hpp:1399
static constexpr Layout layout
Definition views.hpp:1392
VectorViewMut< Scalar > u_box
Definition views.hpp:1404
VEG_INLINE constexpr auto as_const() const noexcept -> QpViewBox< Scalar >
Definition views.hpp:1406
MatrixViewMut< Scalar, layout > A
Definition views.hpp:1397
MatrixViewMut< Scalar, layout > H
Definition views.hpp:1394
VectorViewMut< Scalar > l_box
Definition views.hpp:1403
static constexpr Layout layout
Definition views.hpp:1352
MatrixView< Scalar, layout > H
Definition views.hpp:1354
MatrixView< Scalar, layout > C
Definition views.hpp:1359
MatrixView< Scalar, layout > A
Definition views.hpp:1357
MatrixViewMut< T, layout > A
Definition views.hpp:1375
MatrixViewMut< T, layout > C
Definition views.hpp:1377
VEG_INLINE constexpr auto as_const() const noexcept -> QpView< T >
Definition views.hpp:1380
MatrixViewMut< T, layout > H
Definition views.hpp:1372
static constexpr Layout layout
Definition views.hpp:1370
MatrixView< T, layout > A
Definition views.hpp:1343
MatrixView< T, layout > H
Definition views.hpp:1340
static constexpr Layout layout
Definition views.hpp:1338
MatrixView< T, layout > C
Definition views.hpp:1345
auto operator()(T x) const -> T
Definition views.hpp:1448
auto operator()(Eigen::MatrixBase< D > const &mat) const -> typename D::Scalar
Definition views.hpp:1427
auto operator()(T x, T y) const -> T
Definition views.hpp:1418
auto operator()(T x) const -> T
Definition views.hpp:1439
VEG_INLINE ~Defer() noexcept(noexcept(VEG_FWD(fn)()))
Definition views.hpp:88
Eigen::InnerStride< Eigen::Dynamic > NextColStride
Definition views.hpp:288
static VEG_INLINE auto next_row_stride(isize outer_stride) noexcept -> NextRowStride
Definition views.hpp:289
static VEG_INLINE auto next_col_stride(isize outer_stride) noexcept -> NextColStride
Definition views.hpp:295
static VEG_INLINE constexpr auto offset(T *ptr, isize row, isize col, isize outer_stride) noexcept -> T *
Definition views.hpp:279
static VEG_INLINE void transpose_if_rowmajor(T *ptr, isize dim, isize outer_stride)
Definition views.hpp:303
static VEG_INLINE auto next_row_stride(isize outer_stride) noexcept -> NextRowStride
Definition views.hpp:331
static VEG_INLINE void transpose_if_rowmajor(T *ptr, isize dim, isize outer_stride)
Definition views.hpp:339
static VEG_INLINE constexpr auto offset(T *ptr, isize row, isize col, isize outer_stride) noexcept -> T *
Definition views.hpp:315
Eigen::InnerStride< Eigen::Dynamic > NextRowStride
Definition views.hpp:324
static VEG_INLINE auto next_col_stride(isize outer_stride) noexcept -> NextColStride
Definition views.hpp:325
NoCopy(NoCopy const &)=delete
auto operator=(NoCopy &&) -> NoCopy &=delete
auto operator=(NoCopy const &) -> NoCopy &=delete
proxsuite::linalg::veg::ith< I, Args... > Arg
Definition views.hpp:32
static void fn(T *dest, usize n)
Definition views.hpp:164
VEG_INLINE constexpr auto operator()(Fn fn) const -> Defer< Fn >
Definition views.hpp:95
VEG_INLINE constexpr auto operator()(T const &a, T const &b) const -> T const &
Definition views.hpp:103
VEG_INLINE auto operator()(std::initializer_list< T > list) const -> T
Definition views.hpp:141
VEG_INLINE constexpr auto operator()(T a, T b) const -> T
Definition views.hpp:111