Homec4science
Diffusion lammps (master) Browse

lammps/src/CNEURONET/weigh_container.hh
f9a4fb669e70master

Last Change
Raw File

weigh_container.hh
View Options

#include <valarray>
#include <array>
#include <numeric>
#include <functional>
#include <iterator>
#include <stdexcept>
#include <sstream>
template <size_t Dim, typename T=double>
class DimNumArray {
public:
using index_type = std::array<size_t, Dim>;
DimNumArray(const index_type & shape_={0});
void reshape(const index_type shape_);
inline size_t get_index(const index_type & ijk);
inline T & operator[](const index_type & ijk);
inline T & at(const index_type & ijk);
template <size_t Dir>
class iterator{
public:
inline iterator(size_t start,
DimNumArray & arr, bool end = false);
inline iterator<Dir+1> begin();
inline iterator<Dir+1> end();
inline bool operator!=(const iterator & other);
inline iterator & operator++();
template <size_t ldir=Dir>
inline std::enable_if_t<(ldir<Dim-1), const iterator &>
operator*() const{return this;}
template<size_t ldir=Dir>
inline std::enable_if_t<(ldir==Dim-1), T&> operator*();
inline size_t get_index();
protected:
DimNumArray & arr;
const size_t start;
size_t i;
};
inline iterator<0> begin();
inline iterator<0> end();
void write(FILE *fp) const;
void read (FILE *fp);
private:
index_type shape;
index_type full_strides;
std::valarray<T> values;
};
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
DimNumArray<Dim, T>::DimNumArray(const index_type & shape_):shape(shape_) {
this->reshape(shape_);
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
void DimNumArray<Dim, T>::reshape(const index_type shape_) {
this->shape = shape_;
this->values.resize
(std::accumulate(std::begin(this->shape),
std::end(this->shape),
1,
std::multiplies<size_t>()));
size_t stride = 1;
for (decltype(0) i = Dim-1; i >= 0; --i) {
this->full_strides[i] = stride;
stride *= this->shape[i];
}
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
inline size_t DimNumArray<Dim, T>::get_index(const index_type & ijk){
size_t index{0};
for (size_t i = 0; i < Dim-1; ++i) {
index +=ijk[i];
index *=this->shape[i+1];
}
index += ijk[Dim-1];
return index;
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
inline T & DimNumArray<Dim, T>::operator[](const index_type & ijk) {
return this->values[this->get_index(ijk)];
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
inline T & DimNumArray<Dim, T>::at(const index_type & ijk) {
auto && index = this->get_index(ijk);
auto TODO = index;
bool valid{true};
for (size_t i{0}; i < ijk.size(); ++i) {
valid &= ijk[i]<this->shape[i];
}
if (index < this->values.size() && valid) {
return this->values[index];
} else {
std::stringstream sstr;
auto print = [](const index_type & id) {
std::stringstream tmp;
tmp << "(";
for (size_t i = 0; i<id.size()-1; ++i) {
tmp << id[i] << ", ";
}
tmp << id.back() << ")";
return tmp.str();
};
sstr << "index " << print(ijk)
<< " is out of range for for array of shape "
<< print(this->shape) << ".";
throw std::out_of_range(sstr.str());
}
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
template<size_t Dir>
inline DimNumArray<Dim, T>::iterator<Dir>::iterator(size_t start,
DimNumArray & arr, bool end)
:arr(arr), start(start), i(end ? arr.shape[Dir] : 0) {}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
template<size_t Dir>
size_t
DimNumArray<Dim, T>::iterator<Dir>::get_index(){
return this->start + this->i*this->arr.full_strides[Dir];
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
template<size_t Dir>
inline DimNumArray<Dim, T>::iterator<Dir+1>
DimNumArray<Dim, T>::iterator<Dir>::end(){
return iterator<Dir+1>(this->get_index(),
this->arr, true);
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
template<size_t Dir>
inline DimNumArray<Dim, T>::iterator<Dir+1>
DimNumArray<Dim, T>::iterator<Dir>::begin(){
return iterator<Dir+1>(this->get_index(),
this->arr);
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
template<size_t Dir>
bool
DimNumArray<Dim, T>::iterator<Dir>::operator!=(const iterator & other){
return (this->i != other.i)
|| (this->start != other.start);
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
template<size_t Dir>
inline DimNumArray<Dim, T>::iterator<Dir> &
DimNumArray<Dim, T>::iterator<Dir>::operator++(){
++this->i;
return *this;
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
template<size_t Dir>
template<size_t ldir>
inline std::enable_if_t<(ldir==Dim-1), T&>
DimNumArray<Dim, T>::iterator<Dir>::operator*(){
return this->arr.values[this->get_index()];
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
inline DimNumArray<Dim, T>::iterator<0>
DimNumArray<Dim, T>::begin() {
return iterator<0>(0, *this);
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
inline DimNumArray<Dim, T>::iterator<0>
DimNumArray<Dim, T>::end() {
return iterator<0>(0, *this, true);
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
void DimNumArray<Dim, T>::write(FILE *fp) const {
fwrite(&this->shape[0], sizeof(size_t), Dim, fp);
fwrite(&this->full_strides[0], sizeof(size_t), Dim, fp);
size_t siz = this->values.size();
fwrite(&siz, sizeof(size_t), 1, fp);
fwrite(&this->values[0], sizeof(T), siz, fp);
}
/* ---------------------------------------------------------------------- */
template<size_t Dim, typename T>
void DimNumArray<Dim, T>::read(FILE *fp) {
fread(&this->shape[0], sizeof(size_t), Dim, fp);
fread(&this->full_strides[0], sizeof(size_t), Dim, fp);
size_t siz;
fread(&siz, sizeof(size_t), 1, fp);
fread(&this->values[0], sizeof(T), siz, fp);
}
c4science · Help