Line data Source code
1 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3 : !!!! !!!!
4 : !!!! ParaMonte: Parallel Monte Carlo and Machine Learning Library. !!!!
5 : !!!! !!!!
6 : !!!! Copyright (C) 2012-present, The Computational Data Science Lab !!!!
7 : !!!! !!!!
8 : !!!! This file is part of the ParaMonte library. !!!!
9 : !!!! !!!!
10 : !!!! LICENSE !!!!
11 : !!!! !!!!
12 : !!!! https://github.com/cdslaborg/paramonte/blob/main/LICENSE.md !!!!
13 : !!!! !!!!
14 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
15 : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16 :
17 : !> \brief
18 : !> This file contains procedure implementations of [pm_matrixTrace](@ref pm_matrixTrace).
19 : !>
20 : !> \finmain
21 : !>
22 : !> \author
23 : !> \AmirShahmoradi, September 1, 2017, 12:00 AM, Institute for Computational Engineering and Sciences (ICES), The University of Texas at Austin
24 :
25 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
26 :
27 : ! Define the increment operations for various definitions of trace.
28 : #if getMatTrace_ENABLED
29 : #define TRACE trace
30 : #define OPERATION +
31 : #define GET_LOG(X) X
32 : #elif getMatMulTrace_ENABLED
33 : #define TRACE mulTrace
34 : #define OPERATION *
35 : #define GET_LOG(X) X
36 : #elif getMatMulTraceLog_ENABLED
37 : #define TRACE logMulTrace
38 : #define OPERATION +
39 : #if IK_ENABLED
40 : #define GET_LOG(X) log(real(X, RKD))
41 : #elif CK_ENABLED || RK_ENABLED
42 : #define GET_LOG(X) log(X)
43 : #else
44 : #error "Unrecognized interface."
45 : #endif
46 : #else
47 : #error "Unrecognized interface."
48 : #endif
49 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
50 : #if (getMatTrace_ENABLED || getMatMulTrace_ENABLED || getMatMulTraceLog_ENABLED) && (DEF_ENABLED || RDP_ENABLED) && XXX_ENABLED
51 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
52 :
53 : integer(IK) :: idim
54 306117 : CHECK_ASSERTION(__LINE__, all(0_IK < shape(mat, IK)), SK_"@getMatTrace(): The condition `all(0 < shape(mat, IK))` must hold. shape(mat) = "//getStr(shape(mat, IK)))
55 102039 : CHECK_ASSERTION(__LINE__, size(mat, 1, IK) == size(mat, 2, IK), SK_"@getMatTrace(): The condition `size(mat, 1) == size(mat, 2)` must hold. shape(mat) = "//getStr(shape(mat, IK)))
56 102039 : TRACE = GET_LOG(mat(1, 1))
57 321800 : do idim = 2, size(mat, 1, IK)
58 321800 : TRACE = TRACE OPERATION GET_LOG(mat(idim, idim))
59 : end do
60 :
61 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
62 : #elif (getMatTrace_ENABLED || getMatMulTrace_ENABLED || getMatMulTraceLog_ENABLED) && RFP_ENABLED && (UXD_ENABLED || XLD_ENABLED)
63 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
64 :
65 : integer(IK) :: idim, nrow, ncol
66 0 : nrow = size(mat, 1, IK)
67 0 : ncol = size(mat, 2, IK)
68 0 : CHECK_ASSERTION(__LINE__, all(0_IK < shape(mat, IK)), SK_"@getMatTrace(): The condition `all(0 < shape(mat, IK))` must hold. shape(mat) = "//getStr(shape(mat, IK)))
69 0 : if (ncol < nrow) then ! regular RFP packing, no transposition.
70 0 : if (nrow < 2 * ncol) then ! odd order.
71 0 : CHECK_ASSERTION(__LINE__, nrow == 2 * ncol - 1, SK_"@getMatTrace(): The condition `size(mat, 1) == 2 * size(mat, 2) - 1)` must hold. shape(mat) = "//getStr(shape(mat, IK)))
72 : #if UXD_ENABLED
73 : #define OFFSET_PLUS(X) X +
74 0 : TRACE = GET_LOG(mat(ncol, 1))
75 0 : do idim = 1, ncol - 1
76 0 : TRACE = TRACE OPERATION GET_LOG(mat(ncol + idim, idim))
77 0 : TRACE = TRACE OPERATION GET_LOG(mat(ncol + idim, idim + 1))
78 : end do
79 : #elif XLD_ENABLED
80 : #define OFFSET_PLUS(X)
81 0 : TRACE = GET_LOG(mat(1, 1))
82 0 : do idim = 2, ncol
83 0 : TRACE = TRACE OPERATION GET_LOG(mat(idim - 1, idim))
84 0 : TRACE = TRACE OPERATION GET_LOG(mat(idim, idim))
85 : end do
86 : #endif
87 : else ! even order.
88 0 : CHECK_ASSERTION(__LINE__, nrow == 2 * ncol + 1, SK_"@getMatTrace(): The condition `size(mat, 1) == 2 * size(mat, 2) + 1)` must hold. shape(mat) = "//getStr(shape(mat, IK)))
89 : idim = 1
90 0 : TRACE = GET_LOG(mat(OFFSET_PLUS(ncol) idim, 1))
91 0 : TRACE = TRACE OPERATION GET_LOG(mat(OFFSET_PLUS(ncol) idim + 1, 1))
92 0 : do idim = 2, ncol
93 0 : TRACE = TRACE OPERATION GET_LOG(mat(OFFSET_PLUS(ncol) idim , idim))
94 0 : TRACE = TRACE OPERATION GET_LOG(mat(OFFSET_PLUS(ncol) idim + 1 , idim))
95 : end do
96 : end if
97 0 : elseif (nrow < ncol) then ! transposed RFP packing.
98 0 : if (ncol < 2 * nrow) then ! odd order.
99 0 : CHECK_ASSERTION(__LINE__, ncol == 2 * nrow - 1, SK_"@getMatTrace(): The condition `size(mat, 2) == 2 * size(mat, 1) - 1)` must hold. shape(mat) = "//getStr(shape(mat, IK)))
100 : #if UXD_ENABLED
101 0 : TRACE = GET_LOG(mat(1, nrow))
102 0 : do idim = 1, nrow - 1
103 0 : TRACE = TRACE OPERATION GET_LOG(mat(idim , nrow + idim))
104 0 : TRACE = TRACE OPERATION GET_LOG(mat(idim + 1, nrow + idim))
105 : end do
106 : #elif XLD_ENABLED
107 0 : TRACE = GET_LOG(mat(1, 1))
108 0 : do idim = 2, nrow
109 0 : TRACE = TRACE OPERATION GET_LOG(mat(idim, idim - 1))
110 0 : TRACE = TRACE OPERATION GET_LOG(mat(idim, idim))
111 : end do
112 : #endif
113 : else ! even order.
114 0 : CHECK_ASSERTION(__LINE__, ncol == 2 * nrow + 1, SK_"@getMatTrace(): The condition `size(mat, 2) == 2 * size(mat, 1) + 1)` must hold. shape(mat) = "//getStr(shape(mat, IK)))
115 : idim = 1
116 0 : TRACE = GET_LOG(mat(1, OFFSET_PLUS(nrow) idim))
117 0 : TRACE = TRACE OPERATION GET_LOG(mat(1, OFFSET_PLUS(nrow) idim + 1))
118 0 : do idim = 2, nrow
119 0 : TRACE = TRACE OPERATION GET_LOG(mat(idim, OFFSET_PLUS(nrow) idim ))
120 0 : TRACE = TRACE OPERATION GET_LOG(mat(idim, OFFSET_PLUS(nrow) idim + 1))
121 : end do
122 : end if
123 : else ! nrow == ncol == 1
124 0 : CHECK_ASSERTION(__LINE__, all(1_IK == shape(mat, IK)), SK_"@getMatTrace(): The input RFP `mat` must conform to an upper-triangular subset in RDP packing. shape(mat) = "//getStr(shape(mat, IK)))
125 0 : TRACE = GET_LOG(mat(1, 1))
126 : end if
127 :
128 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
129 : #elif (getMatTrace_ENABLED || getMatMulTrace_ENABLED || getMatMulTraceLog_ENABLED) && LFP_ENABLED && (UXD_ENABLED || XLD_ENABLED)
130 : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
131 :
132 : integer(IK) :: idim, index, ndim, nell
133 15 : nell = size(mat, 1, IK)
134 15 : CHECK_ASSERTION(__LINE__, 0_IK < nell, SK_"@getMatTrace(): The condition `0 < size(mat)` must hold. size(mat) = "//getStr(nell))
135 : !ndim = (getSqrt(nell * 8 + 1) - 1) / 2
136 15 : ndim = (getSqrt(nell * 8 + 1) - 1) / 2
137 60 : CHECK_ASSERTION(__LINE__, nell == ndim * (ndim + 1) / 2, SK_"@getMatTrace(): The input LFP `mat` must conform to an upper-triangular subset in RDP packing. ndim, shape(mat) = "//getStr([ndim, shape(mat, IK)]))
138 15 : TRACE = GET_LOG(mat(1))
139 15 : if (1_IK < nell) then
140 : #if UXD_ENABLED
141 : index = 1
142 0 : do idim = 2, ndim
143 0 : index = index + idim
144 0 : TRACE = TRACE OPERATION GET_LOG(mat(index))
145 : end do
146 : #elif XLD_ENABLED
147 : index = 0
148 58 : do idim = 0, ndim - 2
149 47 : index = index + ndim - idim
150 58 : TRACE = TRACE OPERATION GET_LOG(mat(index))
151 : end do
152 : #endif
153 : end if
154 :
155 : #else
156 : !%%%%%%%%%%%%%%%%%%%%%%%%
157 : #error "Unrecognized interface."
158 : !%%%%%%%%%%%%%%%%%%%%%%%%
159 : #endif
160 : #undef OFFSET_PLUS
161 : #undef OPERATION
162 : #undef GET_LOG
163 : #undef TRACE
|
