TBCI Numerical high perf. C++ Library 2.8.0
vector_extra.h
Go to the documentation of this file.
1
5
6/* KG: 1998/08/17 */
7/* $Id: vector_extra.h,v 1.9.2.17 2022/11/03 17:28:12 garloff Exp $ */
8
9#ifndef TBCI_VECTOR_EXTRA_H
10#define TBCI_VECTOR_EXTRA_H
11
12#include "tbci/vector.h"
13#include "tbci/matrix.h"
14
15// avoid -fguiding-decls (here: unnecessary)
16#if !defined(NO_GD) && !defined(AUTO_DECL)
17# include "tbci/vector_extra_gd.h"
18#endif
19
20NAMESPACE_TBCI
21
22// elementwise multiplication
23INST(template <typename T> class TVector friend TVector<T> emul (const Vector<T>&, const Vector<T>&);)
24template <typename T>
25inline TVector<T> emul (const Vector<T>& a, const Vector<T>& b)
26{
27 const unsigned long dim = b.size();
28 TVector<T> s (dim);
29 BCHK(a.size() != dim, VecErr, product with diff. size Vectors, dim, s);
30 do_vec_vec_mul<T> (dim, &s.setval(0), &a.getcref(0), &b.getcref(0));
31 return s;
32}
33
34INST(template <typename T> class TVector friend TVector<T> emul (TVector<T>, const Vector<T>&);)
35template <typename T>
36inline TVector<T> emul (TVector<T> a, const Vector<T>& b)
37{
38 const unsigned long dim = b.size();
39 BCHK(a.size() != dim, VecErr, product with diff. size Vectors, dim, a);
40 do_vec_mul_vec<T> (dim, &a.setval(0), &b.getcref(0));
41 return a;
42}
43
44INST(template <typename T> class TVector friend TVector<T> emul (const Vector<T>&, TVector<T>);)
45template <typename T>
46inline TVector<T> emul (const Vector<T>& a, TVector<T> b)
47{
48 const unsigned long dim = b.size();
49 BCHK(a.size() != dim, VecErr, product with diff. size Vectors, dim, b);
50 do_vec_mul_vec<T> (dim, &b.setval(0), &a.getcref(0));
51 return b;
52}
53
54INST(template <typename T> class TVector friend TVector<T> emul (TVector<T>, TVector<T>);)
55template <typename T>
56inline TVector<T> emul (TVector<T> a, TVector<T> b)
57{
58 const unsigned long dim = b.size();
59 BCHK(a.size() != dim, VecErr, product with diff. size Vectors, dim, a);
60 do_vec_mul_vec<T> (dim, &a.setval(0), &b.getcref(0));
61 b.destroy (); return a;
62}
63
64INST(template <typename T> class TVector friend TVector<T> cemul (const Vector<T>&, const Vector<T>&);)
65template <typename T>
66inline TVector<T> cemul (const Vector<T>& a, const Vector<T>& b)
67{
68 const unsigned long dim = b.size();
69 TVector<T> s(dim);
70 BCHK(a.size() != dim, VecErr, product with diff. size Vectors, dim, s);
71 do_vec_vec_cmul<T> (dim, &s.setval(0), &a.getcref(0), &b.getcref(0));
72 return s;
73}
74
75INST(template <typename T> class TVector friend TVector<T> cemul (TVector<T>, const Vector<T>&);)
76template <typename T>
77inline TVector<T> cemul (TVector<T> a, const Vector<T>& b)
78{
79 const unsigned long dim = b.size();
80 BCHK(a.size() != dim, VecErr, product with diff. size Vectors, dim, a);
81 do_vec_cmul_vec<T> (dim, &a.setval(0), &b.getcref(0));
82 return a;
83}
84
85INST(template <typename T> class TVector friend TVector<T> cemul (const Vector<T>&, TVector<T>);)
86template <typename T>
87inline TVector<T> cemul (const Vector<T>& a, TVector<T> b)
88{
89 const unsigned long dim = b.size();
90 BCHK(a.size() != dim, VecErr, product with diff. size Vectors, dim, b);
91 do_vec_cmul_vec_inv<T> (dim, &b.setval(0), &a.getcref(0));
92 return b;
93}
94
95INST(template <typename T> class TVector friend TVector<T> cemul (TVector<T>, TVector<T>);)
96template <typename T>
97inline TVector<T> cemul (TVector<T> a, TVector<T> b)
98{
99 const unsigned long dim = b.size();
100 BCHK(a.size() != dim, VecErr, product with diff. size Vectors, dim, a);
101 do_vec_cmul_vec<T> (dim, &a.setval(0), &b.getcref(0));
102 b.destroy (); return a;
103}
104
105
106// elementwise division
107
108INST(template <typename T> class TVector friend TVector<T> ediv (const Vector<T>&, const Vector<T>&);)
109template <typename T>
110inline TVector<T> ediv (const Vector<T>& a, const Vector<T>& b)
111{
112 TVector<T> s (b.size());
113 BCHK(a.size() != b.size(), VecErr, division with diff. size Vectors, b.size(), s);
114 do_vec_vec_div<T> (b.size(), &s.setval(0), &a.getcref(0), &b.getcref(0));
115 return s;
116}
117
118INST(template <typename T> class TVector friend TVector<T> ediv (TVector<T>, const Vector<T>&);)
119template <typename T>
120inline TVector<T> ediv (TVector<T> a, const Vector<T>& b)
121{
122 BCHK(a.size() != b.size(), VecErr, division with diff. size Vectors, b.size(), a);
123 do_vec_div_vec<T> (b.size(), &a.setval(0), &b.getcref(0));
124 return a;
125}
126
127INST(template <typename T> class TVector friend TVector<T> ediv (const Vector<T>&, TVector<T>);)
128template <typename T>
129inline TVector<T> ediv (const Vector<T>& a, TVector<T> b)
130{
131 BCHK(a.size() != b.size(), VecErr, division with diff. size Vectors, b.size(), b);
132 do_vec_div_vec_inv<T> (b.size(), &b.setval(0), &a.getcref(0));
133 return b;
134}
135
136INST(template <typename T> class TVector friend TVector<T> ediv (TVector<T>, TVector<T>);)
137template <typename T>
138inline TVector<T> ediv (TVector<T> a, TVector<T> b)
139{
140 BCHK(a.size() != b.size(), VecErr, division with diff. size Vectors, b.size(), a);
141 do_vec_div_vec<T> (b.size(), &a.setval(0), &b.getcref(0));
142 b.destroy (); return a;
143}
144
145INST(template <typename T> class TVector friend TVector<T> cediv (const Vector<T>&, const Vector<T>&);)
146template <typename T>
147inline TVector<T> cediv (const Vector<T>& a, const Vector<T>& b)
148{
149 TVector<T> s (b.size());
150 BCHK(a.size() != b.size(), VecErr, division with diff. size Vectors, b.size(), s);
151 do_vec_vec_cdiv<T> (b.size(), &s.setval(0), &a.getcref(0), &b.getcref(0));
152 return s;
153}
154
155INST(template <typename T> class TVector friend TVector<T> cediv (TVector<T>, const Vector<T>&);)
156template <typename T>
157inline TVector<T> cediv (TVector<T> a, const Vector<T>& b)
158{
159 BCHK(a.size() != b.size(), VecErr, division with diff. size Vectors, b.size(), a);
160 do_vec_cdiv_vec<T> (b.size(), &a.setval(0), &b.getcref(0));
161 return a;
162}
163
164INST(template <typename T> class TVector friend TVector<T> cediv (const Vector<T>&, TVector<T>);)
165template <typename T>
166inline TVector<T> cediv (const Vector<T>& a, TVector<T> b)
167{
168 BCHK(a.size() != b.size(), VecErr, division with diff. size Vectors, b.size(), b);
169 do_vec_cdiv_vec_inv<T> (a.size(), &b.setval(0), &a.getcref(0));
170 return b;
171}
172
173INST(template <typename T> class TVector friend TVector<T> cediv (TVector<T>, TVector<T>);)
174template <typename T>
175inline TVector<T> cediv (TVector<T> a, TVector<T> b)
176{
177 BCHK(a.size() != b.size(), VecErr, division with diff. size Vectors, b.size(), a);
178 do_vec_cdiv_vec<T> (b.size(), &a.setval(0), &b.getcref(0));
179 b.destroy (); return a;
180}
181
182
183// Now for the combinations with TSVector<T>
184// Not yet optimized with the unrolled kernels from vec_kern_unr_pref.h
185
186INST(template <typename T> class TSVector friend TVector<T> emul (TVector<T>, const TSVector<T>&);)
187template <typename T>
188TVector<T> emul (TVector<T> tv, const TSVector<T>& ts)
189{
190 BCHK(tv.size() != ts.size(), VecErr, Wrong dim in emul, ts.size(), tv);
191 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
192 tv.setval(i) *= ts.get(i);
193 ts.destroy (); return tv;
194}
195
196INST(template <typename T> class TSVector friend TVector<T> emul (const TSVector<T>&, TVector<T>);)
197template <typename T>
198TVector<T> emul (const TSVector<T>& ts, TVector<T> tv)
199{
200 BCHK(tv.size() != ts.size(), VecErr, Wrong dim in emul, tv.size(), tv);
201 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
202 tv.set (ts.get(i) * tv.get(i), i);
203 ts.destroy (); return tv;
204}
205
206
207INST(template <typename T> class TSVector friend TVector<T> emul (const Vector<T>&, const TSVector<T>&);)
208template <typename T>
209TVector<T> emul (const Vector<T>& v, const TSVector<T>& ts)
210{
211 TVector<T> tv;
212 if (ts.mut) { tv.setptr(ts.vecptr()); tv.setsize(ts.size()); }
213 else tv.resize (ts.size());
214 BCHK(v.size() != ts.size(), VecErr, Wrong dim in emul, ts.size(), tv);
215 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
216 tv.set (v.get(i) * ts.get(i), i);
217 return tv;
218}
219
220INST(template <typename T> class TSVector friend TVector<T> emul (const TSVector<T>&, const Vector<T>&);)
221template <typename T>
222TVector<T> emul (const TSVector<T>& ts, const Vector<T>& v)
223{
224 TVector<T> tv;
225 if (ts.mut) { tv.setptr(ts.vecptr()); tv.setsize(ts.size()); }
226 else tv.resize (ts.size());
227 BCHK(v.size() != ts.size(), VecErr, Wrong dim in emul, v.size(), tv);
228 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
229 tv.set (ts.get(i) * v.get(i), i);
230 return tv;
231}
232
233INST(template <typename T> class TSVector friend TVector<T> emul (const TSVector<T>&, const TSVector<T>&);)
234template <typename T>
235TVector<T> emul (const TSVector<T>& ts1, const TSVector<T>& ts2)
236{
237 TVector<T> tv;
238 if (ts1.mut) { tv.setptr(ts1.vecptr()); tv.setsize(ts1.size()); }
239 else if (ts2.mut) { tv.setptr(ts2.vecptr()); tv.setsize(ts2.size()); }
240 else tv.resize (ts1.size());
241 BCHK(ts1.size() != ts2.size(), VecErr, Wrong dim in emul, ts2.size(), tv);
242 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
243 tv.set (ts1.get(i), i);
244 if (ts1.mut && ts2.mut) ts2.destroy ();
245 return tv;
246}
247
248
249INST(template <typename T> class TSVector friend TVector<T> cemul (TVector<T>, const TSVector<T>&);)
250template <typename T>
251TVector<T> cemul (TVector<T> tv, const TSVector<T>& ts)
252{
253 BCHK(tv.size() != ts.size(), VecErr, Wrong dim in emul, ts.size(), tv);
254 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
255 tv.set (CPLX__ conj(tv.get(i)) * ts.get(i), i);
256 ts.destroy (); return tv;
257}
258
259INST(template <typename T> class TSVector friend TVector<T> cemul (const TSVector<T>&, TVector<T>);)
260template <typename T>
261TVector<T> cemul (const TSVector<T>& ts, TVector<T> tv)
262{
263 BCHK(tv.size() != ts.size(), VecErr, Wrong dim in emul, tv.size(), tv);
264 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
265 tv.set (CPLX__ conj (ts.get(i)) * tv.get(i), i);
266 ts.destroy (); return tv;
267}
268
269INST(template <typename T> class TSVector friend TVector<T> cemul (const Vector<T>&, const TSVector<T>&);)
270template <typename T>
271TVector<T> cemul (const Vector<T>& v, const TSVector<T>& ts)
272{
273 TVector<T> tv;
274 if (ts.mut) { tv.setptr(ts.vecptr()); tv.setsize(ts.size()); }
275 else tv.resize (ts.size());
276 BCHK(v.size() != ts.size(), VecErr, Wrong dim in emul, ts.size(), tv);
277 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
278 tv.set (CPLX__ conj(v.get(i)) * ts.get(i), i);
279 return tv;
280}
281
282
283INST(template <typename T> class TSVector friend TVector<T> cemul (const TSVector<T>&, const Vector<T>&);)
284template <typename T>
285TVector<T> cemul (const TSVector<T>& ts, const Vector<T>& v)
286{
287 TVector<T> tv;
288 if (ts.mut) { tv.setptr(ts.vecptr()); tv.setsize(ts.size()); }
289 else tv.resize (ts.size());
290 BCHK(v.size() != ts.size(), VecErr, Wrong dim in emul, v.size(), tv);
291 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
292 tv.set (CPLX__ conj (ts.get(i)) * v.get(i), i);
293 return tv;
294}
295
296INST(template <typename T> class TSVector friend TVector<T> cemul (const TSVector<T>&, const TSVector<T>&);)
297template <typename T>
298TVector<T> cemul (const TSVector<T>& ts1, const TSVector<T>& ts2)
299{
300 TVector<T> tv;
301 if (ts1.mut) { tv.setptr(ts1.vecptr()); tv.setsize(ts1.size()); }
302 else if (ts2.mut) { tv.setptr(ts2.vecptr()); tv.setsize(ts2.size()); }
303 else tv.resize (ts1.size());
304 BCHK(ts1.size() != ts2.size(), VecErr, Wrong dim in emul, ts2.size(), tv);
305 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
306 tv.set (CPLX__ conj (ts1.get(i)), i);
307 if (ts1.mut && ts2.mut) ts2.destroy ();
308 return tv;
309}
310
311
312INST(template <typename T> class TSVector friend TVector<T> ediv (TVector<T>, const TSVector<T>&);)
313template <typename T>
314TVector<T> ediv (TVector<T> tv, const TSVector<T>& ts)
315{
316 BCHK(tv.size() != ts.size(), VecErr, Wrong dim in ediv, ts.size(), tv);
317 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
318 tv.setval(i) /= ts.get(i);
319 ts.destroy (); return tv;
320}
321
322INST(template <typename T> class TSVector friend TVector<T> ediv (const Vector<T>&, const TSVector<T>&);)
323template <typename T>
324TVector<T> ediv (const Vector<T>& v, const TSVector<T>& ts)
325{
326 TVector<T> tv;
327 if (ts.mut) { tv.setptr(ts.vecptr()); tv.setsize(ts.size()); }
328 else tv.resize (ts.size());
329 BCHK(v.size() != ts.size(), VecErr, Wrong dim in ediv, ts.size(), tv);
330 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
331 tv.set (v.get(i) / ts.get(i), i);
332 return tv;
333}
334
335INST(template <typename T> class TSVector friend TVector<T> ediv (const TSVector<T>&, const TSVector<T>&);)
336template <typename T>
337TVector<T> ediv (const TSVector<T>& ts1, const TSVector<T>& ts2)
338{
339 TVector<T> tv;
340 if (ts1.mut) { tv.setptr(ts1.vecptr()); tv.setsize(ts1.size()); }
341 else if (ts2.mut) { tv.setptr(ts2.vecptr()); tv.setsize(ts2.size()); }
342 else tv.resize (ts1.size());
343 BCHK(ts1.size() != ts2.size(), VecErr, Wrong dim in ediv, ts2.size(), tv);
344 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
345 tv.set (ts1.get(i), i);
346 if (ts1.mut && ts2.mut) ts2.destroy ();
347 return tv;
348}
349
350INST(template <typename T> class TSVector friend TVector<T> ediv (const TSVector<T>&, const Vector<T>&);)
351template <typename T>
352TVector<T> ediv (const TSVector<T>& ts, const Vector<T>& v)
353{
354 TVector<T> tv;
355 if (ts.mut) { tv.setptr(ts.vecptr()); tv.setsize(ts.size()); }
356 else tv.resize (ts.size());
357 BCHK(v.size() != ts.size(), VecErr, Wrong dim in ediv, v.size(), tv);
358 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
359 tv.set (ts.get(i) / v.get(i), i);
360 return tv;
361}
362
363INST(template <typename T> class TSVector friend TVector<T> ediv (const TSVector<T>&, TVector<T>);)
364template <typename T>
365TVector<T> ediv (const TSVector<T>& ts, TVector<T> tv)
366{
367 BCHK(tv.size() != ts.size(), VecErr, Wrong dim in ediv, tv.size(), tv);
368 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
369 tv.set (ts.get(i) / tv.get(i), i);
370 ts.destroy (); return tv;
371}
372
373
374INST(template <typename T> class TSVector friend TVector<T> cediv (TVector<T>, const TSVector<T>&);)
375template <typename T>
376TVector<T> cediv (TVector<T> tv, const TSVector<T>& ts)
377{
378 BCHK(tv.size() != ts.size(), VecErr, Wrong dim in ediv, ts.size(), tv);
379 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
380 tv.set (CPLX__ conj(tv.get(i)) / ts.get(i), i);
381 ts.destroy (); return tv;
382}
383
384INST(template <typename T> class TSVector friend TVector<T> cediv (const TSVector<T>&, TVector<T>);)
385template <typename T>
386TVector<T> cediv (const TSVector<T>& ts, TVector<T> tv)
387{
388 BCHK(tv.size() != ts.size(), VecErr, Wrong dim in ediv, tv.size(), tv);
389 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
390 tv.set (CPLX__ conj (ts.get(i)) / tv.get(i), i);
391 ts.destroy (); return tv;
392}
393
394INST(template <typename T> class TSVector friend TVector<T> cediv (const Vector<T>&, const TSVector<T>&);)
395template <typename T>
396TVector<T> cediv (const Vector<T>& v, const TSVector<T>& ts)
397{
398 TVector<T> tv;
399 if (ts.mut) { tv.setptr(ts.vecptr()); tv.setsize(ts.size()); }
400 else tv.resize (ts.size());
401 BCHK(v.size() != ts.size(), VecErr, Wrong dim in ediv, ts.size(), tv);
402 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
403 tv.set (CPLX__ conj(v.get(i)) / ts.get(i), i);
404 return tv;
405}
406
407
408INST(template <typename T> class TSVector friend TVector<T> cediv (const TSVector<T>&, const Vector<T>&);)
409template <typename T>
410TVector<T> cediv (const TSVector<T>& ts, const Vector<T>& v)
411{
412 TVector<T> tv;
413 if (ts.mut) { tv.setptr(ts.vecptr()); tv.setsize(ts.size()); }
414 else tv.resize (ts.size());
415 BCHK(v.size() != ts.size(), VecErr, Wrong dim in ediv, v.size(), tv);
416 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
417 tv.set (CPLX__ conj (ts.get(i)) / v.get(i), i);
418 return tv;
419}
420
421INST(template <typename T> class TSVector friend TVector<T> cediv (const TSVector<T>&, const TSVector<T>&);)
422template <typename T>
423TVector<T> cediv (const TSVector<T>& ts1, const TSVector<T>& ts2)
424{
425 TVector<T> tv;
426 if (ts1.mut) { tv.setptr(ts1.vecptr()); tv.setsize(ts1.size()); }
427 else if (ts2.mut) { tv.setptr(ts2.vecptr()); tv.setsize(ts2.size()); }
428 else tv.resize (ts1.size());
429 BCHK(ts1.size() != ts2.size(), VecErr, Wrong dim in ediv, ts2.size(), tv);
430 for (REGISTER unsigned long i = 0; i < tv.size(); i++)
431 tv.set (CPLX__ conj (ts1.get(i)), i);
432 if (ts1.mut && ts2.mut) ts2.destroy ();
433 return tv;
434}
435
436
437
438INST(template <typename T> class TVector friend TVector<T> ediv (const T&, TVector<T>);)
439template <typename T>
440inline TVector<T> ediv (const T& a, TVector<T> b)
441{
442 STD_SMP_TEMPLATE2C(val_div_vec, b.size(), &b.setval(0), a);
443 return b;
444}
445
446INST(template <typename T> class TVector friend TVector<T> ediv (const T&, const Vector<T>&);)
447template <typename T>
448inline TVector<T> ediv (const T& a, const Vector<T>& b)
449{
450 const unsigned d = b.size();
451 TVector<T> res (d);
452 STD_SMP_TEMPLATE3VC(val_vec_div, d, &res.setval(0), &b.getcref(0), a);
453 return res;
454}
455
456INST(template <typename T> class TSVector friend TVector<T> ediv (const T&, const TSVector<T>&);)
457template <typename T>
458inline TVector<T> ediv (const T& a, const TSVector<T>& ts)
459{
460 TSVector<T> tm (ts); tm.detach ();
461 STD_SMP_TEMPLATE4C(val_svc_div, tm.size(), tm.vecptr(), &ts.getcref(0), a, ts.getfac());
462 tm.getfacref() = (T)1; return TVector<T> (tm);
463}
464
465/* Summing up a vector is seldomly relevant.
466 * TODO: Figure a way to use the algorithms below for other reductions,
467 * such as dot products, fabs, Mat-Vec and Mat-Mat products.
468 * Also create SIMD and parallel versions of it.
469 */
470
471/* http://en.wikipedia.org/wiki/Kahan_summation_algorithm */
472INST(template <typename T> class Vector friend T kahan_sum(const Vector<T>&);)
473template <typename T> T kahan_sum(const Vector<T>& v)
474{
475 T sum = T(0.0); T comp = T(0.0);
476 const unsigned long sz = v.size();
477 for (REGISTER unsigned long i = 0; i < sz; ++i) {
478 T y = v(i) - comp;
479 T t = sum + y;
480 comp = (t - sum) - y;
481 sum = t;
482 }
483 return sum /* - comp */;
484}
485
486/* Improved algo */
487INST(template <typename T> class Vector friend T exact_sum(const Vector<T>&);)
488template <typename T> T exact_sum(const Vector<T>& v)
489{
490 T sum = T(0.0); T comp = T(0.0);
491 const unsigned long sz = v.size();
492 for (REGISTER unsigned long i = 0; i < sz; ++i) {
493 T y = v(i);
494 T t = sum + y;
495 comp += (t - sum) - y;
496 sum = t;
497 }
498 return sum - comp;
499}
500
501// elementwise multiplication / division (with and without taking complex conj.)
502// for matrices
503
504# define _VEC getvec()
505# define _ENDVEC getendvec()
506# define _DIM size()
507# define _ROW rows()
508# define _COL columns()
509# define _FAC getfac()
510
511INST(template <typename T> class TMatrix friend TMatrix<T> emul (const Matrix<T>&, const Matrix<T>&);)
512template <typename T>
513TMatrix<T> emul (const Matrix<T>& a, const Matrix<T>& b)
514{
515 const unsigned int rw = a.rows(), cl = a.columns();
516 BCHK(rw != b.rows(), MatErr, rows differ in emul, rw, TMatrix<T>(0,0));
517 BCHK(cl != b.columns(), MatErr, columns differ in emul, cl, TMatrix<T>(0,0));
518 TMatrix<T> res (rw, cl);
519 do_vec_vec_mul<T> (a._DIM, res._VEC, a._VEC, b._VEC);
520 return res;
521}
522
523INST(template <typename T> class TMatrix friend TMatrix<T> cemul (const Matrix<T>&, const Matrix<T>&);)
524template <typename T>
525TMatrix<T> cemul (const Matrix<T>& a, const Matrix<T>& b)
526{
527 const unsigned int rw = a .rows(), cl = a.columns();
528 BCHK(rw != b.rows(), MatErr, rows differ in cemul, rw, TMatrix<T>(0,0));
529 BCHK(cl != b.columns(), MatErr, columns differ in cemul, cl, TMatrix<T>(0,0));
530 TMatrix<T> res (rw, cl);
531 do_vec_vec_cmul<T> (a._DIM, res._VEC, a._VEC, b._VEC);
532 return res;
533}
534
535INST(template <typename T> class TMatrix friend TMatrix<T> ediv (const Matrix<T>&, const Matrix<T>&);)
536template <typename T>
537TMatrix<T> ediv (const Matrix<T>& a, const Matrix<T>& b)
538{
539 const unsigned int rw = a .rows(), cl = a.columns();
540 BCHK(rw != b.rows(), MatErr, rows differ in ediv, rw, TMatrix<T>(0,0));
541 BCHK(cl != b.columns(), MatErr, columns differ in ediv, cl, TMatrix<T>(0,0));
542 TMatrix<T> res (rw, cl);
543 do_vec_vec_div<T> (a._DIM, res._VEC, a._VEC, b._VEC);
544 return res;
545}
546
547INST(template <typename T> class TMatrix friend TMatrix<T> cediv (const Matrix<T>&, const Matrix<T>&);)
548template <typename T>
549TMatrix<T> cediv (const Matrix<T>& a, const Matrix<T>& b)
550{
551 const unsigned int rw = a .rows(), cl = a.columns();
552 BCHK(rw != b.rows(), MatErr, rows differ in cediv, rw, TMatrix<T>(0,0));
553 BCHK(cl != b.columns(), MatErr, columns differ in cediv, cl, TMatrix<T>(0,0));
554 TMatrix<T> res (rw, cl);
555 do_vec_vec_cdiv<T> (a._DIM, res._VEC, a._VEC, b._VEC);
556 return res;
557}
558
559
560INST(template <typename T> class TMatrix friend TMatrix<T> ediv (const T&, TMatrix<T>);)
561template <typename T>
562inline TMatrix<T> ediv (const T& a, TMatrix<T> b)
563{
564 STD_SMP_TEMPLATE2C(val_div_vec, b._DIM, b._VEC, a);
565 return b;
566}
567
568INST(template <typename T> class TMatrix friend TMatrix<T> ediv (const T&, const Matrix<T>&);)
569template <typename T>
570inline TMatrix<T> ediv (const T& a, const Matrix<T>& b)
571{
572 const unsigned r = b.rows(), c = b.columns();
573 TMatrix<T> res (r, c);
574 STD_SMP_TEMPLATE3VC(val_vec_div, b._DIM, res._VEC, b._VEC, a);
575 return res;
576}
577
578INST(template <typename T> class TSMatrix friend TMatrix<T> ediv (const T&, const TSMatrix<T>&);)
579template <typename T>
580inline TMatrix<T> ediv (const T& a, const TSMatrix<T>& ts)
581{
582 TSMatrix<T> tm (ts); tm.detach ();
583 STD_SMP_TEMPLATE4C(val_svc_div, tm._DIM, tm._VEC, ts._VEC, a, ts._FAC);
584 tm._FAC = (T)1; return TMatrix<T> (tm);
585}
586
587#undef _DIM
588#undef _VEC
589#undef _ENDVEC
590#undef _ROW
591#undef _COL
592#undef _FAC
593
594NAMESPACE_END
595
596#endif /* TBCI_VECTOR_EXTRA_H */
i
int i
Definition LM_fit.h:71
y
doublereal y
Definition TOMS_707.C:27
BCHK
#define BCHK(cond, exc, txt, ind, rtval)
Definition basics.h:575
NAMESPACE_END
#define NAMESPACE_END
Definition basics.h:323
INST
#define INST(x)
Definition basics.h:238
NAMESPACE_TBCI
#define NAMESPACE_TBCI
Definition basics.h:317
REGISTER
#define REGISTER
Definition basics.h:108
T
#define T
Definition bdmatlib.cc:20
BVector::setsize
BVector< T > & setsize(const unsigned long size)
Definition bvector.h:162
BVector::resize
BVector< T > & resize(const BVector< T > &)
Actually it's a resize and copy (some people would expect the assignment op to do this).
Definition bvector.h:361
BVector::setptr
BVector< T > & setptr(T *pointer)
Definition bvector.h:161
MatErr
exception class
Definition matrix.h:54
TSMatrix::detach
void detach(TMatrix< T > *=0)
Definition matrix.h:1239
TSVector::detach
void detach(const T *=0) const
Definition vector.h:1179
TSVector::destroy
void destroy() const
Definition vector.h:1172
TSVector::get
T get(const unsigned long i) const HOT
Definition vector.h:1072
TSVector::getfacref
T & getfacref() const
Definition vector.h:1074
TSVector::mut
bool mut
Definition vector.h:1037
TSVector::vecptr
T *const & vecptr() const
Definition vector.h:1163
TSVector::size
unsigned long size() const
Definition vector.h:1120
TVector
Temporary Base Class Idiom: Class TVector is used for temporary variables.
Definition vector.h:73
TVector::get
tbci_traits< T >::const_refval_type get(const unsigned long i) const
Definition vector.h:190
TVector::size
unsigned long size() const
Definition vector.h:104
TVector::set
T & set(const T &val, const unsigned long i) const
Definition vector.h:194
TVector::setval
T & setval(const unsigned long i) const
Definition vector.h:192
VecErr
exception class
Definition bvector.h:32
conj
NAMESPACE_END NAMESPACE_CPLX TBCI__ cplx< T > conj(const TBCI__ cplx< T > &c)
Definition cplx.h:663
ts
F_TSMatrix< T > ts
Definition f_matrix.h:1052
c
return c
Definition f_matrix.h:760
b
F_TMatrix< T > b
Definition f_matrix.h:736
cediv
FS_Vector< dims, T > cediv(const FS_Vector< dims, T > &f1, const FS_Vector< dims, T > &f2)
Definition fs_vector.h:537
sum
T sum(const FS_Vector< dims, T > &fv)
Definition fs_vector.h:599
ediv
FS_Vector< dims, T > ediv(const FS_Vector< dims, T > &f1, const FS_Vector< dims, T > &f2)
Definition fs_vector.h:528
emul
FS_Vector< dims, T > emul(const FS_Vector< dims, T > &f1, const FS_Vector< dims, T > &f2)
Definition fs_vector.h:510
cemul
FS_Vector< dims, T > cemul(const FS_Vector< dims, T > &f1, const FS_Vector< dims, T > &f2)
Definition fs_vector.h:519
a
const unsigned TMatrix< T > const Matrix< T > * a
Definition matrix_kernels.h:26
res
const unsigned TMatrix< T > * res
Definition matrix_kernels.h:26
CPLX__
STD_SMP_TEMPLATE3VC
#define STD_SMP_TEMPLATE3VC(oper, dm, a1, a2, a3)
Definition vector.h:650
STD_SMP_TEMPLATE2C
#define STD_SMP_TEMPLATE2C(oper, dm, a1, a2)
Definition vector.h:646
STD_SMP_TEMPLATE4C
#define STD_SMP_TEMPLATE4C(oper, dm, a1, a2, a3, a4)
Definition vector.h:656
exact_sum
T exact_sum(const Vector< T > &v)
Definition vector_extra.h:488
kahan_sum
T kahan_sum(const Vector< T > &v)
Definition vector_extra.h:473
Generated by doxygen 1.16.1