PLaSK library

Loading...
Searching...
No Matches
data.hpp
Go to the documentation of this file.
1/*
2 * This file is part of PLaSK (https://plask.app) by Photonics Group at TUL
3 * Copyright (c) 2022 Lodz University of Technology
4 *
5 * This program is free software: you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation, version 3.
8 *
9 * This program is distributed in the hope that it will be useful,
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 * GNU General Public License for more details.
13 */
14#ifndef PLASK__DATA_H
15#define PLASK__DATA_H
16
21#include <iterator>
22#include <algorithm>
23#include <iostream>
24#include <initializer_list>
25#include <atomic>
26#include <type_traits>
27#include <memory> // std::unique_ptr
28#include <cassert>
29#include <type_traits> // std::is_trivially_destructible, std::false_type, std::true_type
30
31#include "memalloc.hpp"
32#include "math.hpp"
33#include "exceptions.hpp"
34#include "log/log.hpp"
35
36#include <boost/type_traits.hpp>
37
38 namespace plask {
39
40 namespace detail {
41
42 template <class T>
43 inline void _construct_array(T* first, T* last, const std::false_type&) {
44 while(first != last) {
45 new(first) T();
46 ++first;
47 }
48 }
49
50 template <class T>
51 inline void _construct_array(T*, T*, const std::true_type&) {
52 }
53
54 template <class T>
55 inline void _destroy_array(T* first, T* last, const std::false_type&) {
56 while(last != first) {
57 --last;
58 last->~T();
59 }
60 }
61
62 template <class T>
63 inline void _destroy_array(T*, T*, const std::true_type&) {
64 }
65
66 template <class T>
67 inline void construct_array(T* first, T* last) {
68 _construct_array(first, last, std::is_trivially_default_constructible<T>());
69 }
70
71 template <class T>
72 inline void destroy_array(T* first, T* last) {
73 _destroy_array(first, last, std::is_trivially_destructible<T>());
74 }
75
76 // We consider dcomplex as POD (do not value-initialize memory)
77
78 template <>
79 inline void construct_array(dcomplex* first, dcomplex* last) {
80 _construct_array(first, last, std::true_type());
81 }
82
83 template <>
84 inline void destroy_array(dcomplex* first, dcomplex* last) {
85 _destroy_array(first, last, std::true_type());
86 }
87
88
90 struct DataVectorGC {
91 // Count is atomic so many threads can increment and decrement it at same time.
92 // If it is 0, it means that there has been only one DataVector object, so probably one thread uses it.
93 std::atomic<unsigned> count;
94
95 std::function<void(void*)>* deleter;
96
97 explicit DataVectorGC(unsigned initial): count(initial), deleter(nullptr) {}
98
99 explicit DataVectorGC(unsigned initial, const std::function<void(void*)>& deleter):
100 count(initial), deleter(new std::function<void(void*)>(deleter)) {}
101
102 explicit DataVectorGC(unsigned initial, std::function<void(void*)>&& deleter):
103 count(initial), deleter(new std::function<void(void*)>(std::forward<std::function<void(void*)>>(deleter))) {}
104
105 void free(void* data) {
106 if (deleter) (*deleter)(data);
107 else aligned_free(data);
108 }
109
110 ~DataVectorGC() { delete deleter; }
111 };
112
113} // namespace detail
114
125 template <typename T>
126 struct DataVector {
127
128 typedef typename std::remove_const<T>::type VT;
129 typedef const T CT;
130
131 typedef detail::DataVectorGC Gc;
132 private:
133
134 std::size_t size_;
135 Gc* gc_;
136 T* data_;
137
139 void dec_ref() { // see http://www.boost.org/doc/libs/1_53_0/doc/html/atomic/usage_examples.html "Reference counting" for optimal memory access description
140 if (gc_ && gc_->count.fetch_sub(1, std::memory_order_release) == 1) {
141 std::atomic_thread_fence(std::memory_order_acquire);
142 detail::destroy_array(data_, data_ + size_);
143 gc_->free(reinterpret_cast<void*>(const_cast<VT*>(data_)));
144 delete gc_;
145 }
146 }
147
149 void inc_ref() {
150 if (gc_) gc_->count.fetch_add(1, std::memory_order_relaxed);
151 }
152
153 friend struct DataVector<VT>;
154 friend struct DataVector<CT>;
155
156 public:
157
158 typedef T value_type;
159
160 typedef T* iterator;
161 typedef const T* const_iterator;
162
164 DataVector(): size_(0), gc_(nullptr), data_(nullptr) {}
165
172 explicit DataVector(std::size_t size): size_(size), gc_(new Gc(1)), data_(aligned_malloc<T>(size)) {
173 detail::construct_array(data_, data_ + size);
174 }
175
181 DataVector(std::size_t size, const T& value): size_(size) {
182 std::unique_ptr<typename std::remove_const<T>::type[], aligned_deleter<T>>
183 data_non_const(aligned_malloc<VT>(size));
184 std::fill_n(data_non_const.get(), size, value); // this may throw, but no memory leak than
185 gc_ = new Gc(1);
186 data_ = data_non_const.release();
187 }
188
193 DataVector(const DataVector<CT>& src): size_(src.size_), gc_(src.gc_), data_(const_cast<T*>(src.data_)) { inc_ref(); }
194
199 DataVector(const DataVector<VT>& src): size_(src.size_), gc_(src.gc_), data_(const_cast<T*>(src.data_)) { inc_ref(); }
200
205 template <typename TS>
206 DataVector(const DataVector<TS>& src): size_(src.size()), gc_(new Gc(1)), data_(aligned_malloc<T>(src.size())) {
207 write_debug("copying data in constructor");
208 std::copy(src.begin(), src.end(), const_cast<VT*>(data_));
209 }
210
216 DataVector<T>& operator=(const DataVector<CT>& M) {
217 const_cast<DataVector<CT>&>(M).inc_ref(); // must be called before dec_ref in case of self-asigment with 1 reference
218 this->dec_ref(); // release old content, this can delete the old data
219 size_ = M.size_;
220 data_ = const_cast<T*>(M.data_);
221 gc_ = M.gc_;
222 return *this;
223 }
224
230 DataVector<T>& operator=(const DataVector<VT>& M) {
231 const_cast<DataVector<VT>&>(M).inc_ref(); // must be called before dec_ref in case of self-asigment with 1 reference
232 this->dec_ref(); // release old content, this can delete the old data
233 size_ = M.size_;
234 data_ = const_cast<T*>(M.data_);
235 gc_ = M.gc_;
236 return *this;
237 }
238
243 DataVector(DataVector<CT>&& src) noexcept: size_(src.size_), gc_(src.gc_), data_(const_cast<T*>(src.data_)) {
244 src.gc_ = nullptr;
245 }
246
251 DataVector(DataVector<VT>&& src) noexcept: size_(src.size_), gc_(src.gc_), data_(const_cast<T*>(src.data_)) {
252 src.gc_ = nullptr;
253 }
254
260 DataVector<T>& operator=(DataVector<T>&& src) noexcept {
261 swap(src);
262 return *this;
263 }
264
270 template <typename TS>
271 DataVector(TS* existing_data, std::size_t size):
272 size_(size), gc_(nullptr), data_(existing_data) {}
273
280 template <typename TS>
281 DataVector(TS* existing_data, std::size_t size, const std::function<void(void*)>& deleter):
282 size_(size), gc_(new Gc(1, deleter)), data_(existing_data) {
283 }
284
291 template <typename TS>
292 DataVector(TS* existing_data, std::size_t size, std::function<void(void*)>&& deleter):
293 size_(size), gc_(new Gc(1, std::forward<std::function<void(void*)>>(deleter))), data_(existing_data) {
294 }
295
300 DataVector(std::initializer_list<T> init): size_(init.size()), gc_(new Gc(1)), data_(aligned_malloc<T>(size_)) {
301 std::copy(init.begin(), init.end(), const_cast<VT*>(data_));
302 }
303
308 template <typename TS>
309 DataVector(std::initializer_list<TS> init): size_(init.size()), gc_(new Gc(1)), data_(aligned_malloc<T>(size_)) {
310 std::copy(init.begin(), init.end(), const_cast<VT*>(data_));
311 }
312
314 ~DataVector() { dec_ref(); }
315
321 void reset() {
322 dec_ref();
323 size_ = 0;
324 gc_ = nullptr;
325 data_ = nullptr;
326 }
327
333 template <typename TS>
334 void reset(TS* existing_data, std::size_t size) {
335 dec_ref();
336 gc_ = nullptr;
337 size_ = size;
338 data_ = existing_data;
339 }
340
347 template <typename TS>
348 void reset(TS* existing_data, std::size_t size, const std::function<void(void*)>& deleter) {
349 dec_ref();
350 gc_ = new Gc(1, deleter);
351 size_ = size;
352 data_ = existing_data;
353 }
354
361 template <typename TS>
362 void reset(TS* existing_data, std::size_t size, std::function<void(void*)>&& deleter) {
363 dec_ref();
364 gc_ = new Gc(1, std::forward<std::function<void(void*)>>(deleter));
365 size_ = size;
366 data_ = existing_data;
367 }
368
377 void reset(std::size_t size) {
378 //TODO to consider: (when clang will be fixed, now it has no std::is_default_constructible but only non-standard std::has_trivial_default_constructor)
379 //if (std::is_default_constructible<T>::value && //this is known at compile time and I believe that compiler optimize-out whole if when it is false
380 // size == size_ && gc_ && gc_->count == 1 && ! gc_->deleter) return;
381 dec_ref();
382 data_ = aligned_malloc<T>(size);
383 detail::construct_array(data_, data_ + size);
384 gc_ = new Gc(1);
385 size_ = size;
386 }
387
395 void reset(std::size_t size, const T& value) {
396 std::unique_ptr<VT[], aligned_deleter<T>>
397 data_non_const(aligned_malloc<VT>(size));
398 std::fill_n(data_non_const.get(), size, value); // this may throw, than our data will not be changed
399 dec_ref();
400 gc_ = new Gc(1); //this also may throw
401 data_ = data_non_const.release();
402 size_ = size;
403 }
404
411 template <typename InIterT>
412 void reset(InIterT begin, InIterT end) {
413 std::unique_ptr<VT[], aligned_deleter<T>>
414 data_non_const(aligned_malloc<VT>(std::distance(begin,end)));
415 std::copy(begin, end, data_non_const.get()); // this may throw, and than our vector will not be changed
416 dec_ref();
417 gc_ = new Gc(1);
418 size_ = std::distance(begin, end);
419 data_ = data_non_const.release();
420 }
421
426 const_iterator begin() const { return data_; }
427
432 iterator begin() { return data_; }
433
438 const_iterator end() const { return data_ + size_; }
439
444 iterator end() { return data_ + size_; }
445
447 std::size_t size() const { return size_; }
448
450 const T* data() const { return data_; }
451
453 T* data() { return data_; }
454
460 const T& operator[](std::size_t n) const { assert(n < size_); return data_[n]; }
461
467 T& operator[](std::size_t n) { assert(n < size_); return data_[n]; }
468
470 operator bool() const { return data_ != nullptr; }
471
476 DataVector<VT> copy() const {
477 DataVector<VT> new_data(size_);
478 std::copy(begin(), end(), new_data.begin());
479 return new_data;
480 }
481
486 bool unique() const {
487 return (gc_ != nullptr) && (gc_->count == 1);
488 }
489
494 DataVector<VT> remove_const() const {
495 DataVector<VT> result(const_cast<VT*>(this->data()), this->size());
496 result.gc_ = this->gc_;
497 result.inc_ref();
498 return result;
499 }
500
505 DataVector<VT> claim() const {
506 return (unique() && !gc_->deleter)? remove_const() : copy();
507 }
508
513 void swap(DataVector<T>& other) noexcept {
514 std::swap(size_, other.size_);
515 std::swap(gc_, other.gc_);
516 std::swap(data_, other.data_);
517 }
518
523 template <typename O>
524 void fill(const O& value) {
525 std::fill_n(data_, size_, value);
526 }
527
534 inline DataVector<T> getSubarrayRef(std::size_t begin_index, std::size_t subarray_size) {
535 assert(begin_index + subarray_size <= size_);
536 return DataVector<T>(data_ + begin_index, subarray_size);
537 }
538
545 inline DataVector<const T> getSubarrayRef(std::size_t begin_index, std::size_t subarray_size) const {
546 assert(begin_index + subarray_size <= size_);
547 return DataVector<const T>(data_ + begin_index, subarray_size);
548 }
549
556 inline DataVector<T> getSubarrayCopy(std::size_t begin_index, std::size_t subarray_size) const {
557 return getSubarrayRef(begin_index, subarray_size).copy();
558 }
559
565 inline DataVector<T> getSubrangeRef(iterator begin, iterator end) {
566 assert(this->begin() <= begin && begin <= end && end <= this->end());
567 return DataVector<T>(begin, end - begin);
568 }
569
575 inline DataVector<const T> getSubrangeRef(const_iterator begin, const_iterator end) const {
576 assert(this->begin() <= begin && begin <= end && end <= this->end());
577 return DataVector<const T>(begin, end-begin);
578 }
579
585 inline DataVector<T> getSubrangeCopy(const_iterator begin, const_iterator end) {
586 assert(this->begin() <= begin && begin <= end && end <= this->end());
587 return getSubrangeRef(begin, end).copy();
588 }
589};
590
596 template<class T1, class T2> inline
597 bool operator== ( DataVector<T1> const& a, DataVector<T2> const& b)
598{ return a.size() == b.size() && std::equal(a.begin(), a.end(), b.begin()); }
599
605 template<class T1, class T2> inline
606 bool operator!= ( DataVector<T1> const& a, DataVector<T2> const& b) { return !(a==b); }
607
613 template<class T1, class T2> inline
614 bool operator< ( DataVector<T1> const& a, DataVector<T2> const& b)
615{ return std::lexicographical_compare(a.begin(), a.end(), b.begin(), b.end()); }
616
622 template<class T1, class T2> inline
623 bool operator> ( DataVector<T1> const& a, DataVector<T2> const& b) { return b < a; }
624
630 template<class T1, class T2> inline
631 bool operator<= ( DataVector<T1> const& a, DataVector<T2> const& b) { return !(b < a); }
632
638 template<class T1, class T2> inline
639 bool operator>= ( DataVector<T1> const& a, DataVector<T2> const& b) { return !(a < b); }
640
647 template<class T>
648std::ostream& operator<<(std::ostream& out, DataVector<T> const& to_print) {
649 out << '[';
650 return print_seq(out, to_print.begin(), to_print.end()) << ']';
651}
652
659 template<class T, class S>
660 DataVector<T>& operator+=(DataVector<T>& to_inc, DataVector<S> const& inc_val) {
661 if (to_inc.size() != inc_val.size())
662 throw DataError("data vectors sizes differ ([{0}] += [%2])", to_inc.size(), inc_val.size());
663 for (std::size_t i = 0; i < to_inc.size(); ++i)
664 to_inc[i] += inc_val[i];
665 return to_inc;
666}
667
674 template<class T, class S>
675 DataVector<T>& operator-=(DataVector<T>& to_dec, DataVector<S> const& dec_val) {
676 if (to_dec.size() != dec_val.size())
677 throw DataError("data vectors sizes differ ([{0}] -= [%2])", to_dec.size(), dec_val.size());
678 for (std::size_t i = 0; i < to_dec.size(); ++i)
679 to_dec[i] -= dec_val[i];
680 return to_dec;
681}
682
689 template<typename T, typename S>
690 DataVector<T>& operator*=(DataVector<T>& vec, S a) {
691 for (std::size_t i = 0; i < vec.size(); ++i)
692 vec[i] *= a;
693 return vec;
694}
695
702 template<typename T, typename S>
703 DataVector<T>& operator/=(DataVector<T>& vec, S a) {
704 auto ia = 1. / a;
705 for (std::size_t i = 0; i < vec.size(); ++i)
706 vec[i] *= ia;
707 return vec;
708}
709
716 template<typename T1, typename T2>
717 DataVector<typename std::remove_cv<decltype(T1()+T2())>::type> operator+(DataVector<T1> const& vec1, DataVector<T2> const& vec2) {
718 if (vec1.size() != vec2.size())
719 throw DataError("data vectors sizes differ ([{0}] + [%2])", vec1.size(), vec2.size());
720 DataVector<typename std::remove_cv<decltype(T1()+T2())>::type> result(vec1.size());
721 for (std::size_t i = 0; i < vec1.size(); ++i)
722 result[i] = vec1[i] + vec2[i];
723 return result;
724}
725
732 template<typename T1, typename T2>
733 DataVector<typename std::remove_cv<decltype(T1()-T2())>::type> operator-(DataVector<T1> const& vec1, DataVector<T2> const& vec2) {
734 if (vec1.size() != vec2.size())
735 throw DataError("data vectors sizes differ ([{0}] - [%2])", vec1.size(), vec2.size());
736 DataVector<typename std::remove_cv<decltype(T1()-T2())>::type> result(vec1.size());
737 for (std::size_t i = 0; i < vec1.size(); ++i)
738 result [i] = vec1[i] - vec2[i];
739 return result;
740}
741
747 template<class T>
748 DataVector<typename std::remove_cv<T>::type> operator-(DataVector<T> const& vec) {
749 DataVector<typename std::remove_cv<T>::type> result(vec.size());
750 for (std::size_t i = 0; i < vec.size(); ++i)
751 result [i] = -vec[i];
752 return result;
753}
754
761 template<typename T, typename S>
762 DataVector<typename std::remove_cv<decltype(T()*S())>::type> operator*(DataVector<T> const& vec, S a) {
763 DataVector<typename std::remove_cv<decltype(T()*S())>::type> result(vec.size());
764 for (std::size_t i = 0; i < vec.size(); ++i)
765 result[i] = vec[i] * a;
766 return result;
767}
768
775 template<typename T, typename S>
776 DataVector<typename std::remove_cv<decltype(S()*T())>::type> operator*(S a, DataVector<T> const& vec) {
777 DataVector<typename std::remove_cv<decltype(S()*T())>::type> result(vec.size());
778 for (std::size_t i = 0; i < vec.size(); ++i)
779 result[i] = a * vec[i];
780 return result;
781}
782
789 template<typename T, typename S>
790 DataVector<typename std::remove_cv<decltype(T()*(1./S()))>::type> operator/(DataVector<T> const& vec, S a) {
791 auto ia = 1. / a;
792 DataVector<typename std::remove_cv<decltype(T()*(1./S()))>::type> result(vec.size());
793 for (std::size_t i = 0; i < vec.size(); ++i)
794 result[i] = vec[i] * ia;
795 return result;
796}
797
798
799
806 template <class T>
807 typename std::remove_cv<T>::type accumulate(const DataVector<T>& to_accum, typename std::remove_cv<T>::type initial=typename std::remove_cv<T>::type()) {
808 for (std::size_t i = 0; i < to_accum.size(); ++i)
809 initial += to_accum[i];
810 return initial;
811}
812
817 template <class T>
818T average(const DataVector<T>& v) {
819 return accumulate(v) / v.size();
820}
821
828 template<typename RT, typename T>
829 inline DataVector<RT> const_data_cast(const DataVector<T>& src) {
830 return src.remove_const();
831}
832
833/*
834 PLASK_API_EXTERN_TEMPLATE_STRUCT(DataVector<double>)
835 PLASK_API_EXTERN_TEMPLATE_STRUCT(DataVector<const double>)
836 PLASK_API_EXTERN_TEMPLATE_STRUCT(PLASK_API DataVector<std::complex<double>>)
837 PLASK_API_EXTERN_TEMPLATE_STRUCT(DataVector<const std::complex<double>>)
838*/
839
840} // namespace plask
841
842 namespace std {
843 template <typename T>
844 void swap(plask::DataVector<T>& s1, plask::DataVector<T>& s2) noexcept {
845 s1.swap(s2);
846 }
847} // namespace std
848
849#if FMT_VERSION >= 90000
850 template <typename T> struct fmt::formatter<plask::DataVector<T>> : ostream_formatter {};
851#endif
852
853#endif // PLASK__DATA_H