ergo
matrix_proxy.h
Go to the documentation of this file.
1/* Ergo, version 3.8.2, a program for linear scaling electronic structure
2 * calculations.
3 * Copyright (C) 2023 Elias Rudberg, Emanuel H. Rubensson, Pawel Salek,
4 * and Anastasia Kruchinina.
5 *
6 * This program is free software: you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation, either version 3 of the License, or
9 * (at your option) any later version.
10 *
11 * This program is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program. If not, see <http://www.gnu.org/licenses/>.
18 *
19 * Primary academic reference:
20 * Ergo: An open-source program for linear-scaling electronic structure
21 * calculations,
22 * Elias Rudberg, Emanuel H. Rubensson, Pawel Salek, and Anastasia
23 * Kruchinina,
24 * SoftwareX 7, 107 (2018),
25 * <http://dx.doi.org/10.1016/j.softx.2018.03.005>
26 *
27 * For further information about Ergo, see <http://www.ergoscf.org>.
28 */
29
41#ifndef MAT_MATRIX_PROXY
42#define MAT_MATRIX_PROXY
43
44namespace mat {
45 /*********** New code */
50 template<typename TX, typename TY>
51 struct XY {
52 TX const & A;
53 TY const & B;
54 bool const tA;
55 bool const tB;
56 XY(TX const & AA, TY const & BB,
57 bool const tAA = false, bool const tBB = false)
58 :A(AA), B(BB), tA(tAA), tB(tBB)
59 {}
60 };
61
66 template<typename TX, typename TY, typename TZ>
67 struct XYZ {
68 TX const & A;
69 TY const & B;
70 TZ const & C;
71 bool const tA;
72 bool const tB;
73 bool const tC;
74 XYZ(TX const & AA, TY const & BB, TZ const & CC,
75 bool const tAA = false,
76 bool const tBB = false,
77 bool const tCC = false)
78 :A(AA), B(BB), C(CC), tA(tAA), tB(tBB), tC(tCC)
79 {}
80 };
81
87 template<typename TX, typename TY, typename TZ, typename TU, typename TV>
88 struct XYZpUV {
89 TX const & A;
90 TY const & B;
91 TZ const & C;
92 TU const & D;
93 TV const & E;
94 bool const tA;
95 bool const tB;
96 bool const tC;
97 bool const tD;
98 bool const tE;
99 XYZpUV(TX const & AA, TY const & BB, TZ const & CC,
100 TU const & DD, TV const & EE,
101 bool const tAA = false,
102 bool const tBB = false,
103 bool const tCC = false,
104 bool const tDD = false,
105 bool const tEE = false)
106 :A(AA), B(BB), C(CC), D(DD), E(EE),
107 tA(tAA), tB(tBB), tC(tCC), tD(tDD), tE(tEE)
108 {}
109 };
110
111
117 template<typename TX>
118 struct Xtrans {
119 TX const & A;
120 bool const tA;
121 explicit Xtrans(TX const & AA, bool const tAA = false)
122 :A(AA), tA(tAA)
123 {}
124 };
125
130 template<typename TX>
131 inline Xtrans<TX> transpose(TX const & A) {
132 return Xtrans<TX>(A,true);
133 }
141 template<typename TX>
143 return Xtrans<TX>(xtrans.A, !xtrans.tA);
144 }
145
146 /* Some operators */
156 template<typename TX, typename TY>
158 Xtrans<TY> const & trBB) {
159 return XY<TX, TY>(trAA.A, trBB.A, trAA.tA, trBB.tA);
160 }
161
171 template<typename TX, typename TY>
172 inline XY<TX, TY> operator*(TX const & AA,
173 Xtrans<TY> const & trBB) {
174 return XY<TX, TY>(AA, trBB.A, false, trBB.tA);
175 }
176
186 template<typename TX, typename TY>
188 TY const & BB) {
189 return XY<TX, TY>(trAA.A, BB, trAA.tA, false);
190 }
191
199 template<typename TX, typename TY>
200 inline XY<TX, TY> operator*(TX const & AA,
201 TY const & BB) {
202 return XY<TX, TY>(AA, BB, false, false);
203 }
204
212 template<typename TX, typename TY, typename TZ>
213 inline XYZ<TX, TY, TZ>
215 return XYZ<TX, TY, TZ>(AB.A, AB.B, trCC.A, AB.tA, AB.tB, trCC.tA);
216 }
217
223 template<typename TX, typename TY, typename TZ>
224 inline XYZ<TX, TY, TZ>
225 operator*(XY<TX, TY> const & AB, TZ const & CC) {
226 return XYZ<TX, TY, TZ>(AB.A, AB.B, CC, AB.tA, AB.tB, false);
227 }
228
235 template<typename TX, typename TY, typename TZ, typename TU, typename TV>
238 return XYZpUV<TX, TY, TZ, TU, TV>(ABC.A, ABC.B, ABC.C, DE.A, DE.B, ABC.tA, ABC.tB, ABC.tC, DE.tA, DE.tB);
239 }
240
245 template<typename TX, typename TY>
246 struct XpY {
247 const TX& A;
248 const TY& B;
249 XpY(const TX& AA,const TY& BB)
250 :A(AA),B(BB)
251 {}
252 };
256 template<typename TX, typename TY>
257 inline XpY<TX, TY> operator+(TX const & AA, TY const & BB) {
258 return XpY<TX, TY>(AA, BB);
259 }
260
265 template<typename TX, typename TY>
266 struct XmY {
267 const TX& A;
268 const TY& B;
269 XmY(const TX& AA,const TY& BB)
270 :A(AA),B(BB)
271 {}
272 };
276 template<typename TX, typename TY>
277 inline XmY<TX, TY> operator-(TX const & AA, TY const & BB) {
278 return XmY<TX, TY>(AA, BB);
279 }
280
281
282 /************* New code ends */
283
284
285#if 0
286 template<class MAT>
287 struct M2 {
288 const MAT& A;
289 M2(const MAT& AA)
290 :A(AA)
291 {}
292 };
293
294 template<class MAT>
295 inline M2<MAT> square(const MAT& A) {
296 return M2<MAT>(A);
297 }
298
299 template<class SCAL, class MAT>
300 struct SM2 {
301 const SCAL alpha;
302 const MAT& A;
303 SM2(const MAT& AA, const SCAL a = 1)
304 : A(AA), alpha(a)
305 {}
306 SM2(const M2<MAT>& m2)
307 :A(m2.A), alpha(1)
308 {}
309 };
310
311 template<class SCAL, class MAT>
312 inline SM2<SCAL, MAT>
313 operator*(const SCAL s, const M2<MAT>& m2) {
314 return SM2<SCAL, MAT>(m2.A, s);
315 }
316
317
318
319
320 template<class MAT>
321 struct MT {
322 const MAT& A;
323 const bool tA;
324 MT(const MAT& AA, const bool tAA = false)
325 :A(AA), tA(tAA)
326 {}
327 };
328
329 template<class MAT>
330 inline MT<MAT> transpose(const MAT& A) {
331 return MT<MAT>(A,true);
332 }
333 template<class MAT>
334 inline MT<MAT> transpose(const MT<MAT>& mt) {
335 return MT<MAT>(mt.A, !mt.tA);
336 }
337
338
339 template<class SCAL, class MAT>
340 struct SM {
341 const SCAL alpha;
342 const MAT& A;
343 const bool tA;
344 SM(const MAT& AA, const SCAL scalar = 1, const bool tAA = false)
345 :A(AA),alpha(scalar), tA(tAA)
346 {}
347 };
348
349 template<class SCAL, class MAT>
350 inline SM<SCAL, MAT> operator*(const SCAL scalar, const MT<MAT>& mta) {
351 return SM<SCAL, MAT>(mta.A,scalar, mta.tA);
352 }
353
354 template<class SCAL, class MAT>
355 inline SM<SCAL, MAT> operator*(const SCAL scalar, const MAT& AA) {
356 return SM<SCAL, MAT>(AA, scalar, false);
357 }
358
359
360
361 template<class MAT, class MATB = MAT>
362 struct MM {
363 const MAT& A;
364 const MATB& B;
365 const bool tA;
366 const bool tB;
367 MM(const MAT& AA,const MATB& BB, const bool tAA, const bool tBB)
368 :A(AA),B(BB), tA(tAA), tB(tBB)
369 {}
370 };
371
372 template<class MAT, class MATB = MAT>
373 struct MpM {
374 const MAT& A;
375 const MATB& B;
376 MpM(const MAT& AA,const MATB& BB)
377 :A(AA),B(BB)
378 {}
379 };
380
381 template<class MAT, class MATB>
382 inline MpM<MAT, MATB> operator+(const MAT& AA, const MATB& BB) {
383 return MpM<MAT, MATB>(AA, BB);
384 }
385
386 /*
387 template<class MAT, class MATB>
388 inline MM<MAT, MATB> operator*(const MT<MAT>& mta, const MT<MATB>& mtb) {
389 return MM<MAT, MATB>(mta.A, mtb.A, mta.tA, mtb.tA);
390 }
391 */
392 /*
393 template<class MAT, class MATB>
394 inline MM<MAT, MATB> operator*(const MAT& AA, const MT<MATB>& mtb) {
395 return MM<MAT, MATB>(AA, mtb.A, false, mtb.tA);
396 }
397 template<class MAT, class MATB>
398 inline MM<MAT, MATB> operator*(const MT<MAT>& mta, const MATB& BB) {
399 return MM<MAT, MATB>(mta.A, BB, mta.tA, false);
400 }
401 template<class MAT, class MATB>
402 inline MM<MAT, MATB> operator*(const MAT& AA, const MATB& BB) {
403 return MM<MAT, MATB>(AA, BB, false, false);
404 }
405 */
406
407 template<class SCAL, class MAT, class MATB = MAT>
408 struct SMM {
409 const SCAL alpha;
410 const MAT& A;
411 const MATB& B;
412 const bool tA;
413 const bool tB;
414 SMM(const MM<MAT, MATB>& mm)
415 :A(mm.A),B(mm.B),alpha(1), tA(mm.tA), tB(mm.tB)
416 {}
417 SMM(const MAT& AA,const MATB& BB,
418 const bool tAA, const bool tBB,
419 const SCAL a = 1)
420 :A(AA), B(BB), tA(tAA), tB(tBB), alpha(a)
421 {}
422 };
423
424 template<class SCAL, class MAT, class MATB>
426 operator*(const SM<SCAL, MAT>& sm,const MT<MATB>& mtb) {
427 return SMM<SCAL, MAT, MATB>(sm.A, mtb.A, sm.tA, mtb.tA, sm.alpha);
428 }
429
430 template<class SCAL, class MAT, class MATB>
432 operator*(const SM<SCAL, MAT>& sm,const MATB& BB) {
433 return SMM<SCAL, MAT, MATB>(sm.A, BB, sm.tA, false, sm.alpha);
434 }
435
436 template<class SCAL, class MATC, class MATA = MATC, class MATB = MATC>
437 struct SMMpSM {
438 const SCAL alpha;
439 const MATA& A;
440 const MATB& B;
441 const SCAL beta;
442 const MATC& C;
443 const bool tA;
444 const bool tB;
445 SMMpSM(const MATA& AA, const MATB& BB, const MATC& CC,
446 const bool tAA, const bool tBB,
447 const SCAL a=1, const SCAL b=1)
448 :A(AA), B(BB), C(CC), alpha(a), beta(b), tA(tAA), tB(tBB)
449 {}
450 };
451
452 template<class SCAL, class MATC, class MATA, class MATB>
456 (smm.A, smm.B, sm.A, smm.tA, smm.tB, smm.alpha, sm.alpha);
457 }
458#if 0
459
460 template<class SCAL, class MATC, class MATA, class MATB>
464 (smm.A, smm.B, CC, smm.tA, smm.tB, smm.alpha, 1);
465 }
466 template<class SCAL, class MATC, class MATA, class MATB>
470 (mm.A, mm.B, sm.A, mm.tA, mm.tB, 1, sm.alpha);
471 }
472#endif
473
474
475 template<class SCAL, class MAT>
476 struct SM2pSM {
477 const SCAL alpha;
478 const MAT& A;
479 const SCAL beta;
480 const MAT& C;
481 SM2pSM(const MAT& AA, const MAT& CC, const SCAL a = 1, const SCAL b = 0)
482 : A(AA), alpha(a), C(CC), beta(b)
483 {}
484 };
485
486 template<class SCAL, class MAT>
487 inline SM2pSM<SCAL, MAT>
489 return SM2pSM<SCAL, MAT>(sm2.A, sm.A, sm2.alpha, sm.alpha);
490 }
491
492
493
494 /* Done so far with new transpose */
495
496template<class MAT>
497 struct MMpM {
498 const MAT& A;
499 const MAT& B;
500 const MAT& C;
501 MMpM(const MAT& AA, const MAT& BB, const MAT& CC)
502 :A(AA),B(BB),C(CC)
503 {}
504 };
505
506
507 template<class SCAL, class MAT>
508 struct SMpSM {
509 const SCAL alpha, beta;
510 const MAT& A, B;
511 SMpSM(const MAT& AA, const MAT& BB,
512 const SCAL scalar_a=1, const SCAL scalar_b=1)
513 :A(AA), B(BB), alpha(scalar_a), beta(scalar_b)
514 {}
515 };
516
517 template<class SCAL, class MAT>
518 inline SMpSM<SCAL, MAT>
520 return SMpSM<SCAL, MAT>(sm1.A, sm2.A, sm1.alpha, sm2.alpha);
521 }
522
523 /*
524 template<class MAT>
525 struct MpM {
526 const MAT& A;
527 const MAT& B;
528 MpM(const MAT& AA,const MAT& BB)
529 :A(AA),B(BB)
530 {}
531 };
532
533 template<class MAT>
534 inline MpM<MAT> operator+(const MAT& A, const MAT& B) {
535 return MpM<MAT>(A,B);
536 }
537 */
538 template<class MAT>
539 struct MmM {
540 const MAT& A;
541 const MAT& B;
542 MmM(const MAT& AA,const MAT& BB)
543 :A(AA),B(BB)
544 {}
545 };
546
547 template<class MAT>
548 inline MmM<MAT> operator-(const MAT& A, const MAT& B) {
549 return MmM<MAT>(A,B);
550 }
551
552
553
554
555
556
557 /*onodig finns redan för SMM
558 template<class MAT>
559 inline MMpM<MAT> operator+(const MM<MAT>& mm, const MAT& CC)
560 {
561 return MMpM<MAT>(mm.A,mm.B,CC);
562 }*/
563
564
565 /*Maste ligga i arvda klassen!!*/
566 /*
567 Matrix::Matrix(const sMMmul& mm)
568 :nrofrows(mm.A.nrofrows),nrofcols(mm.B.nrofcols)
569 {
570 this.multiply(mm.A,mm.B,*this,mm.tA,mm.tB,mm.alpha,0);
571 }
572 Matrix::Matrix(const sMMmulsMadd& mm)
573 :nrofrows(mm.A.nrofrows),nrofcols(mm.B.nrofcols)
574 {
575 this->multiply(mm.A,mm.B,mm.C,mm.tA,mm.tB,mm.alpha,mm.beta);
576 }
577
578 */
579#endif
580} /* end namespace mat */
581#endif
#define B
#define A
Definition allocate.cc:39
XmY< TX, TY > operator-(TX const &AA, TY const &BB)
Substraction of two objects of type TX and TY.
Definition matrix_proxy.h:277
XYZpUV< TX, TY, TZ, TU, TV > operator+(XYZ< TX, TY, TZ > const &ABC, XY< TU, TV > const &DE)
Addition of two multiplication proxys XYZ and XY.
Definition matrix_proxy.h:237
static Treal getMachineEpsilon()
Definition matInclude.h:147
Xtrans< TX > transpose(TX const &A)
Transposition.
Definition matrix_proxy.h:131
XY< TX, TY > operator*(Xtrans< TX > const &trAA, Xtrans< TY > const &trBB)
Multiplication of two transposition proxys holding objects of type TX and TY respectively.
Definition matrix_proxy.h:157
This proxy expresses the result of multiplication of three objects, of possibly different types,...
Definition matrix_proxy.h:67
TZ const & C
Definition matrix_proxy.h:70
bool const tB
Definition matrix_proxy.h:72
bool const tC
Definition matrix_proxy.h:73
bool const tA
Definition matrix_proxy.h:71
TY const & B
Definition matrix_proxy.h:69
XYZ(TX const &AA, TY const &BB, TZ const &CC, bool const tAA=false, bool const tBB=false, bool const tCC=false)
Definition matrix_proxy.h:74
TX const & A
Definition matrix_proxy.h:68
This proxy expresses the result of multiplication of three objects added to two other multiplied obje...
Definition matrix_proxy.h:88
bool const tC
Definition matrix_proxy.h:96
bool const tE
Definition matrix_proxy.h:98
TZ const & C
Definition matrix_proxy.h:91
TY const & B
Definition matrix_proxy.h:90
bool const tA
Definition matrix_proxy.h:94
TU const & D
Definition matrix_proxy.h:92
TX const & A
Definition matrix_proxy.h:89
bool const tD
Definition matrix_proxy.h:97
TV const & E
Definition matrix_proxy.h:93
bool const tB
Definition matrix_proxy.h:95
XYZpUV(TX const &AA, TY const &BB, TZ const &CC, TU const &DD, TV const &EE, bool const tAA=false, bool const tBB=false, bool const tCC=false, bool const tDD=false, bool const tEE=false)
Definition matrix_proxy.h:99
This proxy expresses the result of multiplication of two objects, of possibly different types,...
Definition matrix_proxy.h:51
TY const & B
Definition matrix_proxy.h:53
TX const & A
Definition matrix_proxy.h:52
bool const tA
Definition matrix_proxy.h:54
bool const tB
Definition matrix_proxy.h:55
XY(TX const &AA, TY const &BB, bool const tAA=false, bool const tBB=false)
Definition matrix_proxy.h:56
This proxy expresses the result of substraction of two objects, of possibly different types,...
Definition matrix_proxy.h:266
XmY(const TX &AA, const TY &BB)
Definition matrix_proxy.h:269
const TY & B
Definition matrix_proxy.h:268
const TX & A
Definition matrix_proxy.h:267
This proxy expresses the result of addition of two objects, of possibly different types,...
Definition matrix_proxy.h:246
const TX & A
Definition matrix_proxy.h:247
XpY(const TX &AA, const TY &BB)
Definition matrix_proxy.h:249
const TY & B
Definition matrix_proxy.h:248
This proxy expresses the result of transposition of an object of type TX.
Definition matrix_proxy.h:118
TX const & A
Definition matrix_proxy.h:119
Xtrans(TX const &AA, bool const tAA=false)
Definition matrix_proxy.h:121
bool const tA
Definition matrix_proxy.h:120