fml  0.1-0
Fused Matrix Library
parmat.hh
1 // This file is part of fml which is released under the Boost Software
2 // License, Version 1.0. See accompanying file LICENSE or copy at
3 // https://www.boost.org/LICENSE_1_0.txt
4 
5 #ifndef FML_PAR_INTERNALS_PARMAT_H
6 #define FML_PAR_INTERNALS_PARMAT_H
7 #pragma once
8 
9 
10 #include <cmath>
11 #include <cstdlib>
12 #include <cstring>
13 #include <random>
14 
15 #include "../../_internals/rand.hh"
16 #include "../../_internals/types.hh"
17 #include "../comm.hh"
18 
19 
20 namespace fml
21 {
22  template <class MAT, class VEC, typename REAL>
23  class parmat
24  {
25  public:
26  parmat(){};
27  parmat(comm mpi_comm, MAT &data_);
28  parmat(parmat<MAT, VEC, REAL> &&x);
29 
30  void resize(len_global_t nrows, len_t ncols);
31 
32  void print(uint8_t ndigits=4, bool add_final_blank=true);
33  void info();
34 
35  void fill_zero();
36  void fill_val(const REAL v);
37  void fill_linspace(const REAL start, const REAL stop);
38  void fill_eye();
39  void fill_diag(const VEC &d);
40  void fill_runif(const uint32_t seed, const REAL min=0, const REAL max=1);
41  void fill_runif(const REAL min=0, const REAL max=1);
42  void fill_rnorm(const uint32_t seed, const REAL mean=0, const REAL sd=1);
43  void fill_rnorm(const REAL mean=0, const REAL sd=1);
44 
45  // void diag(cpuvec<REAL> &v);
46  // void antidiag(cpuvec<REAL> &v);
47  void scale(const REAL s);
48  // void rev_rows();
49  void rev_cols();
50 
51  bool any_inf() const;
52  bool any_nan() const;
53 
54  REAL get(const len_global_t i) const;
55  REAL get(const len_global_t i, const len_t j) const;
56  void set(const len_global_t i, const REAL v);
57  void set(const len_global_t i, const len_t j, const REAL v);
58  // void get_row(const len_global_t i, VEC &v) const;
59  // void get_col(const len_global_t j, VEC &v) const;
60 
61  bool operator==(const parmat<MAT, VEC, REAL> &x) const;
62  bool operator!=(const parmat<MAT, VEC, REAL> &x) const;
63  parmat<MAT, VEC, REAL>& operator=(const parmat<MAT, VEC, REAL> &x);
64 
65  len_global_t nrows() const {return m_global;};
66  len_local_t nrows_local() const {return data.nrows();};
67  len_local_t ncols() const {return data.ncols();};
68  len_global_t nrows_before() const {return nb4;};
69  comm get_comm() const {return r;};
70  const MAT& data_obj() const {return data;};
71  MAT& data_obj() {return data;};
72 
73  protected:
74  MAT data;
75  len_global_t m_global;
76  comm r;
77  len_global_t nb4;
78 
79  void num_preceding_rows();
80  len_t get_local_dim();
81  void check_index(const len_global_t i) const;
82  void check_index(const len_global_t i, const len_t j) const;
83  };
84 }
85 
86 
87 
88 template <class MAT, class VEC, typename REAL>
89 fml::parmat<MAT, VEC, REAL>::parmat(fml::comm mpi_comm, MAT &data_)
90 {
91  r = mpi_comm;
92  data = data_;
93 
94  m_global = (len_global_t) data.nrows();
95  r.allreduce(1, &(m_global));
96  num_preceding_rows();
97 }
98 
99 
100 
101 template <class MAT, class VEC, typename REAL>
102 fml::parmat<MAT, VEC, REAL>::parmat(fml::parmat<MAT, VEC, REAL> &&x)
103 {
104  this->data = x.data_obj();
105  this->m_global = x.nrows();
106  this->r = x.get_comm();
107  this->nb4 = x.nrows_before();
108 }
109 
110 
111 
112 template <class MAT, class VEC, typename REAL>
113 void fml::parmat<MAT, VEC, REAL>::resize(len_global_t nrows, len_t ncols)
114 {
115  this->m_global = nrows;
116  len_t m_local = this->get_local_dim();
117 
118  this->data.resize(m_local, ncols);
119  num_preceding_rows();
120 }
121 
122 
123 
124 template <class MAT, class VEC, typename REAL>
125 void fml::parmat<MAT, VEC, REAL>::info()
126 {
127  r.printf(0, "# parmat");
128  r.printf(0, " %" PRIu64 "x%d", m_global, data.ncols());
129  r.printf(0, " type=%s", typeid(REAL).name());
130  r.printf(0, "\n");
131 }
132 
133 
134 
135 template <class MAT, class VEC, typename REAL>
136 void fml::parmat<MAT, VEC, REAL>::fill_zero()
137 {
138  data.fill_zero();
139 }
140 
141 
142 
143 template <class MAT, class VEC, typename REAL>
144 void fml::parmat<MAT, VEC, REAL>::fill_val(const REAL v)
145 {
146  data.fill_val(v);
147 }
148 
149 
150 
151 template <class MAT, class VEC, typename REAL>
152 void fml::parmat<MAT, VEC, REAL>::fill_runif(const uint32_t seed, const REAL min, const REAL max)
153 {
154  data.fill_runif(seed, min, max);
155 }
156 
157 template <class MAT, class VEC, typename REAL>
158 void fml::parmat<MAT, VEC, REAL>::fill_runif(const REAL min, const REAL max)
159 {
160  uint32_t seed = fml::rand::get_seed() + r.rank();
161  data.fill_runif(seed, min, max);
162 }
163 
164 
165 
166 template <class MAT, class VEC, typename REAL>
167 void fml::parmat<MAT, VEC, REAL>::fill_rnorm(const uint32_t seed, const REAL mean, const REAL sd)
168 {
169  data.fill_rnorm(seed, mean, sd);
170 }
171 
172 template <class MAT, class VEC, typename REAL>
173 void fml::parmat<MAT, VEC, REAL>::fill_rnorm(const REAL mean, const REAL sd)
174 {
175  uint32_t seed = fml::rand::get_seed() + r.rank();
176  data.fill_rnorm(seed, mean, sd);
177 }
178 
179 
180 
181 template <class MAT, class VEC, typename REAL>
182 void fml::parmat<MAT, VEC, REAL>::scale(const REAL v)
183 {
184  data.scale(v);
185 }
186 
187 
188 
189 template <class MAT, class VEC, typename REAL>
190 void fml::parmat<MAT, VEC, REAL>::rev_cols()
191 {
192  data.rev_cols();
193 }
194 
195 
196 
197 template <class MAT, class VEC, typename REAL>
198 bool fml::parmat<MAT, VEC, REAL>::any_inf() const
199 {
200  int ret = (int) data.any_inf();
201  r.allreduce(1, &ret);
202  return (bool) ret;
203 }
204 
205 
206 
207 template <class MAT, class VEC, typename REAL>
208 bool fml::parmat<MAT, VEC, REAL>::any_nan() const
209 {
210  int ret = (int) data.any_nan();
211  r.allreduce(1, &ret);
212  return (bool) ret;
213 }
214 
215 
216 
217 template <class MAT, class VEC, typename REAL>
218 REAL fml::parmat<MAT, VEC, REAL>::get(const len_global_t i) const
219 {
220  check_index(i);
221 
222  len_global_t row = i % m_global;
223  len_t j = (len_t) (i / m_global);
224 
225  REAL ret;
226  if (row >= nb4 && row < nb4+data.nrows())
227  ret = data.get(row-nb4, j);
228  else
229  ret = 0;
230 
231  r.allreduce(1, &ret);
232  return ret;
233 }
234 
235 template <class MAT, class VEC, typename REAL>
236 REAL fml::parmat<MAT, VEC, REAL>::get(const len_global_t i, const len_t j) const
237 {
238  check_index(i, j);
239 
240  REAL ret;
241  if (i >= nb4 && i < nb4+data.nrows())
242  ret = data.get(i-nb4, j);
243  else
244  ret = 0;
245 
246  r.allreduce(1, &ret);
247  return ret;
248 }
249 
250 template <class MAT, class VEC, typename REAL>
251 void fml::parmat<MAT, VEC, REAL>::set(const len_global_t i, const REAL v)
252 {
253  check_index(i);
254 
255  len_global_t row = i % m_global;
256  len_t j = (len_t) (i / m_global);
257 
258  if (row >= nb4 && row < nb4+data.nrows())
259  data.set(row-nb4, j, v);
260 }
261 
262 template <class MAT, class VEC, typename REAL>
263 void fml::parmat<MAT, VEC, REAL>::set(const len_global_t i, const len_t j, const REAL v)
264 {
265  check_index(i, j);
266 
267  if (i >= nb4 && i < nb4+data.nrows())
268  data.set(i-nb4, j, v);
269 }
270 
271 
272 
273 template <class MAT, class VEC, typename REAL>
274 bool fml::parmat<MAT, VEC, REAL>::operator==(const fml::parmat<MAT, VEC, REAL> &x) const
275 {
276  int neq_count = (int) (data != x.data_obj());
277  r.allreduce(1, &neq_count);
278 
279  return (neq_count == 0);
280 }
281 
282 template <class MAT, class VEC, typename REAL>
283 bool fml::parmat<MAT, VEC, REAL>::operator!=(const fml::parmat<MAT, VEC, REAL> &x) const
284 {
285  return !(*this == x);
286 }
287 
288 
289 
290 template <class MAT, class VEC, typename REAL>
291 fml::parmat<MAT, VEC, REAL>& fml::parmat<MAT, VEC, REAL>::operator=(const fml::parmat<MAT, VEC, REAL> &x)
292 {
293  this->data = x.data_obj();
294  this->m_global = x.nrows();
295  this->r = x.get_comm();
296  this->nb4 = x.nrows_before();
297 }
298 
299 
300 // -----------------------------------------------------------------------------
301 // private
302 // -----------------------------------------------------------------------------
303 
304 template <class MAT, class VEC, typename REAL>
305 void fml::parmat<MAT, VEC, REAL>::num_preceding_rows()
306 {
307  int myrank = r.rank();
308  int size = r.size();;
309 
310  nb4 = 0;
311  len_t m_local = data.nrows();
312 
313  for (int rank=1; rank<size; rank++)
314  {
315  if (myrank == (rank - 1))
316  {
317  len_global_t nb4_send = nb4 + ((len_global_t) m_local);
318  r.send(1, &nb4_send, rank);
319  }
320  else if (myrank == rank)
321  {
322  len_global_t nr_prev_rank;
323  r.recv(1, &nr_prev_rank, rank-1);
324 
325  nb4 += nr_prev_rank;
326  }
327  }
328 }
329 
330 
331 
332 template <class MAT, class VEC, typename REAL>
333 len_t fml::parmat<MAT, VEC, REAL>::get_local_dim()
334 {
335  len_t local = m_global / r.size();
336  len_t rem = (len_t) (m_global - (len_global_t) local*r.size());
337  if (r.rank()+1 <= rem)
338  local++;
339 
340  return local;
341 }
342 
343 
344 
345 template <class MAT, class VEC, typename REAL>
346 void fml::parmat<MAT, VEC, REAL>::check_index(const len_global_t i) const
347 {
348  if (i < 0 || i >= (m_global * data.ncols()))
349  throw std::runtime_error("index out of bounds");
350 }
351 
352 template <class MAT, class VEC, typename REAL>
353 void fml::parmat<MAT, VEC, REAL>::check_index(const len_global_t i, const len_t j) const
354 {
355  if (i < 0 || i >= m_global || j < 0 || j >= data.ncols())
356  throw std::runtime_error("index out of bounds");
357 }
358 
359 
360 #endif
fml::parmat
Definition: parmat.hh:23
fml::comm::allreduce
void allreduce(int n, T *data, MPI_Op op=MPI_SUM) const
Sum reduce operation across all processes in the MPI communicator.
Definition: comm.hh:357
fml::comm
MPI communicator data and helpers.
Definition: comm.hh:24
fml
Core namespace.
Definition: dimops.hh:10