ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_matrixUpdate::setMatUpdate Interface Reference

Return the result of arbitrary-rank Symmetric/Hermitian updates to triangular matrices of type integer, complex, and real of arbitrary type-kind parameters. More...

Detailed Description

Return the result of arbitrary-rank Symmetric/Hermitian updates to triangular matrices of type integer, complex, and real of arbitrary type-kind parameters.

see the documentation of pm_matrixUpdate for more details.

Parameters
[in,out]mat: The input/output contiguous matrix of arbitrary shape (:,:) of the same type and kind as matA, containing the subset of matrix \(\ms{mat}\) in the triangular matrix update.
[in]subset: The input scalar that can be either,
  1. the constant uppDia implying that the upper-diagonal triangular block of mat should be used for storage and updating while the lower triangle remains intact.
  2. the constant lowDia implying that the lower-diagonal triangular block of mat should be used for storage and updating while the upper triangle remains intact.
This argument is merely a convenience to differentiate the different procedure functionalities within this generic interface.
[in]matA: The input contiguous matrix of arbitrary shape (:,:) of,
  1. type integer of kind any supported by the processor (e.g., IK, IK8, IK16, IK32, or IK64), or
  2. type complex of kind any supported by the processor (e.g., CK, CK32, CK64, or CK128), or
  3. type real of kind any supported by the processor (e.g., RK, RK32, RK64, or RK128),
containing the subset of matrix \(\ms{matA}\) in the triangular matrix update.
[in]operationA: The input scalar parameter that can be,
  1. the constant transSymm if matA is of type integer, complex, or real, implying the use of the Symmetric form of the specified subset of matA in the triangular matrix update.
  2. the constant transHerm if matA is of type complex and alpha and beta are of type real, implying the use of the Hermitian transpose form of the specified subset of matA in the triangular matrix update.
Specifying this argument changes the shape of the subset of matA used in the triangular matrix update. See the description of the input argument ndum.
This argument is merely a convenience to differentiate the different procedure functionalities within this generic interface.
(optional. If missing, the specified subset of matA will be used as is, without transposition.)
[in]alpha: The input scalar containing the coefficient \(\alpha\) in the triangular matrix update.
  1. If matA is of type integer or real, then alpha must of the same type and kind as matA.
  2. If matA is of type complex, then alpha must be of the type complex or real of the same kind as matA.
    1. If alpha is of type complex, then the resulting matrix update will be Symmetric of the form,

      \begin{align*} & \ms{mat} \leftarrow \alpha \ms{matA}~~ \ms{matA}^T + \beta \ms{mat} &&\text{if operationA is missing.} \\ & \ms{mat} \leftarrow \alpha \ms{matA}^T \ms{matA}~~ + \beta \ms{mat} &&\text{if operationA = transSymm} \end{align*}

    2. If alpha is of type real, then the resulting matrix update will be Hermitian of the form,

      \begin{align*} & \ms{mat} \leftarrow \alpha \ms{matA}~~ \ms{matA}^H + \beta \ms{mat} &&\text{if operationA is missing.} \\ & \ms{mat} \leftarrow \alpha \ms{matA}^H \ms{matA}~~ + \beta \ms{mat} &&\text{if operationA = transHerm} \end{align*}

[in]beta: The input scalar containing the coefficient \(\beta\) in the triangular matrix update.
  1. If matA is of type integer or real, then beta must of the same type and kind as matA.
  2. If matA is of type complex, then beta must be of the type complex or real of the same kind as matA.
    1. If beta is of type complex, then the resulting matrix update will be Symmetric of the form,

      \begin{align*} & \ms{mat} \leftarrow \beta \ms{matA}~~ \ms{matA}^T + \beta \ms{mat} &&\text{if operationA is missing.} \\ & \ms{mat} \leftarrow \beta \ms{matA}^T \ms{matA}~~ + \beta \ms{mat} &&\text{if operationA = transSymm} \end{align*}

    2. If beta is of type real, then the resulting matrix update will be Hermitian of the form,

      \begin{align*} & \ms{mat} \leftarrow \beta \ms{matA}~~ \ms{matA}^H + \beta \ms{mat} &&\text{if operationA is missing.} \\ & \ms{mat} \leftarrow \beta \ms{matA}^H \ms{matA}~~ + \beta \ms{mat} &&\text{if operationA = transHerm} \end{align*}

[in]ndim: The input non-negative scalar of type integer of default kind IK, containing the number of rows/columns (i.e., the rank) of mat used in the update mat(roff + 1 : roff + ndim, coff + 1 : coff + ndim).
[in]ndum: The input non-negative scalar of type integer of default kind IK, containing the number of dummy rows or columns of matA used in the triangular matrix update.
  1. For matrix matA,
    1. If the input argument operationA is missing, then the subset matA(roffA + 1 : roffA + ndim, coffA + 1 : coffA + ndum) is used in the triangular matrix update.
    2. If the input argument operationA is transSymm, then the subset transpose(matA(roffA + 1 : roffA + ndum, coffA + 1 : coffA + ndim) is used in the triangular matrix update.
    3. If the input argument operationA is transHerm, then the subset conjg(transpose(matA(roffA + 1 : roffA + ndum, coffA + 1 : coffA + ndim)) is used in the triangular matrix update.
[in]roff: The input non-negative scalar of type integer of default kind IK containing the offset from the first row of the input matrix mat, such that mat(1 + roff, 1 + coff) marks the top-left corner of the block of mat used in the triangular matrix update.
[in]coff: The input non-negative scalar of type integer of default kind IK containing the offset from the first column of the input matrix mat, such that mat(1 + roff, 1 + coff) marks the top-left corner of the block of mat used in the triangular matrix update.
[in]roffA: The input non-negative scalar of type integer of default kind IK containing the offset from the first row of the input matrix matA, such that matA(1 + roffA, 1 + coffA) marks the top-left corner of the block of matA used in the triangular matrix update.
[in]coffA: The input non-negative scalar of type integer of default kind IK containing the offset from the first column of the input matrix matA, such that matA(1 + roffA, 1 + coffA) marks the top-left corner of the block of matA used in the triangular matrix update.


Possible calling interfaces

use pm_matrixUpdate, only: transSymm, transHerm ! Possible values of `operationA`
use pm_matrixUpdate, only: uppDia, lowDia ! subset: upper-diagonal, lower-diagonal
call setMatUpdate(mat, subset, matA, operationA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
call setMatUpdate(matA = transSymm, mat, subset, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
call setMatUpdate(matA = transHerm, mat, subset, alpha, beta, ndim, ndum, roff, coff, roffA, coffA) ! complex(*) :: matA, real(kind(matA)) :: alpha, beta
!
Return the result of arbitrary-rank Symmetric/Hermitian updates to triangular matrices of type intege...
This module contains procedures and generic interfaces relevant to arbitrary-rank updates to vectors,...
Warning
The condition alphaim == 0 must hold when class = hermitian.
The condition all(0 <= [roffA, coffA]) must hold for the corresponding input arguments.
The condition all(0 <= [roff, coff]) must hold for the corresponding input arguments.
The condition all(0_IK <= [ndim, ndim, ndum]) must hold for the corresponding input arguments.
The condition all([roffA + ndim, coffA + ndum] <= shape(matA)) must hold for the corresponding input arguments when the input argument operationA is missing.
The condition all([roffA + ndum, coffA + ndim] <= shape(matA)) must hold for the corresponding input arguments when the input argument operationA is present.
This condition is verified only if the library is built with the preprocessor macro CHECK_ENABLED=1.
The pure procedure(s) documented herein become impure when the ParaMonte library is compiled with preprocessor macro CHECK_ENABLED=1.
By default, these procedures are pure in release build and impure in debug and testing builds.
Developer Remark:
The optional input argument operationA in the current interface was originally designed be either missing or set to trans standing for transposed.
While such convention offers a cleaner interface, it was later replaced with the two possible options transSymm and transHerm.
Although the new interface is uglier and verbose, it allows the possibility of extending this interface in the future without breaking the existing interface.
The redundancy also offers an automatic compile-time check against semantic bugs, given that alpha, beta must be of type real when transHerm transposition is requested.
BLAS/LAPACK equivalent:
The procedures under discussion combine, modernize, and extend the interface and functionalities of Version 3.11 of BLAS/LAPACK routine(s): SSYRK, DSYRK, CSYRK, ZSYRK, CHERK, and ZHERK.
Notably, the interfaces are also extended to support matrices of type integer of arbitrary kinds.

BLAS/LAPACK argument correspondence
The following table demonstrates the equivalence of the arguments of the relevant BLAS/LAPACK routines with the arguments of the corresponding interface in the ParaMonte library.

BLAS/LAPACK interface arguments setMatUpdate interface arguments
uplo = "U" subset = uppDia
uplo = "L" subset = lowDia
n ndim
k ndum
trans = "N" operationA argument is missing.
trans = "T" operationA = transSymm
trans = "C" operationA = transHerm
alpha alpha
a = A(i,j) matA = A, roffA = i - 1, coffA = j - 1
lda NONE (passed implicitly).
beta beta
c = C(i,j) mat = C, roff = i - 1, coff = j - 1
ldc NONE (passed implicitly).
See also
uppDia
lowDia
transSymm
transHerm
setMatUpdateR1


Example usage

1program example
2
3 use pm_kind, only: SK, IK, LK
4 use pm_kind, only: RKG => RKS ! all processor type kinds are supported.
5 use pm_kind, only: CKG => CKS ! all processor type kinds are supported.
6 use pm_matrixUpdate, only: lowDia, uppDia
7 use pm_matrixUpdate, only: symmetric, hermitian
8 use pm_matrixUpdate, only: nothing, trans
10 use pm_io, only: getFormat
11 use pm_io, only: display_type
12
13 implicit none
14
15 type(display_type) :: disp
16
17 real(RKG) , parameter :: DUM = -huge(0._RKG)
18 complex(CKG) , parameter :: CMPLX_DUMM = cmplx(-huge(0._CKG), -huge(0._CKG), CKG)
19 integer(IK) :: ndim, ndum, roffC, coffC, roffA, coffA
20 character(:, SK), allocatable :: cform
21 cform = getFormat([cmplx(0., 0., CKG)], ed = SK_'f', signed = .true.)
22
23 disp = display_type(file = "main.out.F90")
24
25 block
26
27 real(RKG) :: alpha, beta
28 real(RKG), allocatable :: RefC(:,:), triC(:,:), matC(:,:), matA(:,:), VecX(:), VecY(:)
29
30 matA = reshape( [ real(RKG) :: 0.0, 8.0 &
31 , 1.0, 9.0 &
32 , 2.0, 10.0 &
33 , 3.0, 11.0 &
34 , 4.0, 12.0 &
35 , 5.0, 13.0 &
36 , 6.0, 14.0 &
37 , 7.0, 15.0 &
38 , DUM, DUM &
39 ], shape = [9, 2], order = [2, 1])
40 triC = reshape( [ real(RKG) :: 0.0, 1.0, 3.0, 6.0, 10.0, 15.0, 21.0, 28.0 &
41 , DUM, 2.0, 4.0, 7.0, 11.0, 16.0, 22.0, 29.0 &
42 , DUM, DUM, 5.0, 8.0, 12.0, 17.0, 23.0, 30.0 &
43 , DUM, DUM, DUM, 9.0, 13.0, 18.0, 24.0, 31.0 &
44 , DUM, DUM, DUM, DUM, 14.0, 19.0, 25.0, 32.0 &
45 , DUM, DUM, DUM, DUM, DUM, 20.0, 26.0, 33.0 &
46 , DUM, DUM, DUM, DUM, DUM, DUM, 27.0, 34.0 &
47 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, 35.0 &
48 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
49 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
50 ], shape = [10, 8], order = [2, 1])
51 RefC = reshape( [ real(RKG):: 64.0, 73.0, 83.0, 94.0, 106.0, 119.0, 133.0, 148.0 &
52 , DUM, 84.0, 96.0, 109.0, 123.0, 138.0, 154.0, 171.0 &
53 , DUM, DUM, 109.0, 124.0, 140.0, 157.0, 175.0, 194.0 &
54 , DUM, DUM, DUM, 139.0, 157.0, 176.0, 196.0, 217.0 &
55 , DUM, DUM, DUM, DUM, 174.0, 195.0, 217.0, 240.0 &
56 , DUM, DUM, DUM, DUM, DUM, 214.0, 238.0, 263.0 &
57 , DUM, DUM, DUM, DUM, DUM, DUM, 259.0, 286.0 &
58 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, 309.0 &
59 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
60 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
61 ], shape = [10, 8], order = [2, 1])
62 matC = triC
63 call disp%skip()
64 call disp%show("matA")
65 call disp%show( matA )
66 call disp%show("matC")
67 call disp%show( matC )
68 call disp%show("alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 2; roffC = 0; coffC = 0; roffA = 0; coffA = 0;")
69 alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 2; roffC = 0; coffC = 0; roffA = 0; coffA = 0;
70 call disp%show("call setMatUpdate(matC, symmetric, uppDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?SYRK contiguous interface.")
71 call setMatUpdate(matC, symmetric, uppDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?SYRK contiguous interface.
72 call disp%show("matC")
73 call disp%show( matC )
74 call disp%show("matC - RefC")
75 call disp%show( matC - RefC )
76 call disp%skip()
77 matC = triC
78 call disp%show("call setMatUpdate(matC(1:8,1:8), symmetric, uppDia, matA(1:8,1:2), nothing) ! BLAS 3 ?SYRK simplified interface.")
79 call setMatUpdate(matC(1:8,1:8), symmetric, uppDia, matA(1:8,1:2), nothing) ! BLAS 3 ?SYRK simplified interface.
80 call disp%show("matC")
81 call disp%show( matC )
82 call disp%show("matC - RefC")
83 call disp%show( matC - RefC )
84 call disp%skip()
85
86 matA = reshape( [ real(RKG) :: 0.0, 3.0, 6.0, 9.0, 12.0, 15.0, 18.0, 21.0 &
87 , 1.0, 4.0, 7.0, 10.0, 13.0, 16.0, 19.0, 22.0 &
88 , 2.0, 5.0, 8.0, 11.0, 14.0, 17.0, 20.0, 23.0 &
89 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
90 ], shape = [4, 8], order = [2, 1])
91 triC = reshape( [ real(RKG) :: 0.0, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
92 , 1.0, 8.0, DUM, DUM, DUM, DUM, DUM, DUM &
93 , 2.0, 9.0, 15.0, DUM, DUM, DUM, DUM, DUM &
94 , 3.0, 10.0, 16.0, 21.0, DUM, DUM, DUM, DUM &
95 , 4.0, 11.0, 17.0, 22.0, 26.0, DUM, DUM, DUM &
96 , 5.0, 12.0, 18.0, 23.0, 27.0, 30.0, DUM, DUM &
97 , 6.0, 13.0, 19.0, 24.0, 28.0, 31.0, 33.0, DUM &
98 , 7.0, 14.0, 20.0, 25.0, 29.0, 32.0, 34.0, 35.0 &
99 ], shape = [8, 8], order = [2, 1])
100 RefC = reshape( [ real(RKG) :: 5.0, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
101 , 15.0, 58.0, DUM, DUM, DUM, DUM, DUM, DUM &
102 , 25.0, 95.0, 164.0, DUM, DUM, DUM, DUM, DUM &
103 , 35.0, 132.0, 228.0, 323.0, DUM, DUM, DUM, DUM &
104 , 45.0, 169.0, 292.0, 414.0, 535.0, DUM, DUM, DUM &
105 , 55.0, 206.0, 356.0, 505.0, 653.0, 800.0, DUM, DUM &
106 , 65.0, 243.0, 420.0, 596.0, 771.0, 945.0, 1118.0, DUM &
107 , 75.0, 280.0, 484.0, 687.0, 889.0, 1090.0, 1290.0, 1489.0 &
108 ], shape = [8, 8], order = [2, 1])
109 matC = triC
110 call disp%skip()
111 call disp%show("matA")
112 call disp%show( matA )
113 call disp%show("matC")
114 call disp%show( matC )
115 call disp%show("alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 3; roffA = 0; coffA = 0; roffC = 0; coffC = 0;")
116 alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 3; roffA = 0; coffA = 0; roffC = 0; coffC = 0;
117 call disp%show("call setMatUpdate(matC, symmetric, lowDia, matA, trans, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?SYRK contiguous interface.")
118 call setMatUpdate(matC, symmetric, lowDia, matA, trans, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?SYRK contiguous interface.
119 call disp%show("matC")
120 call disp%show( matC )
121 call disp%show("matC - RefC")
122 call disp%show( matC - RefC )
123 call disp%skip()
124 matC = triC
125 call disp%show("call setMatUpdate(matC(1:8, 1:8), symmetric, lowDia, matA(1:3, 1:8), trans) ! BLAS 3 ?SYRK simplified interface.")
126 call setMatUpdate(matC(1:8, 1:8), symmetric, lowDia, matA(1:3, 1:8), trans) ! BLAS 3 ?SYRK simplified interface.
127 call disp%show("matC")
128 call disp%show( matC )
129 call disp%show("matC - RefC")
130 call disp%show( matC - RefC )
131 call disp%skip()
132
133 end block
134
135 call disp%skip
136 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
137 call disp%show("!Complex Symmetric matrix update.")
138 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
139 call disp%skip
140
141 block
142
143 complex(CKG) :: alpha, beta
144 complex(CKG), allocatable :: RefC(:,:), triC(:,:), matC(:,:), matA(:,:), VecX(:), VecY(:)
145
146 matA = reshape( [ complex(CKG) :: (2.0, 0.0), (3.0, 2.0), (4.0, 1.0), (1.0, 7.0), (0.0, 0.0) &
147 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0), (2.0, 4.0), (1.0, 2.0) &
148 , (1.0, 3.0), (2.0, 1.0), (6.0, 0.0), (3.0, 2.0), (2.0, 2.0) &
149 ], shape = [3, 5], order = [2, 1])
150 triC = reshape( [ complex(CKG) :: (2.0, 1.0), (1.0, 9.0), (4.0, 5.0) &
151 , CMPLX_DUMM, (3.0, 1.0), (6.0, 7.0) &
152 , CMPLX_DUMM, CMPLX_DUMM, (8.0, 1.0) &
153 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
154 ], shape = [4, 3], order = [2, 1])
155 RefC = reshape( [ complex(CKG) :: (-57.0, 13.0), (-63.0, 79.0), (-24.0, 70.0) &
156 , CMPLX_DUMM, (-28.0, 90.0), (-55.0, 103.0) &
157 , CMPLX_DUMM, CMPLX_DUMM, (13.0, 75.0) &
158 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
159 ], shape = [4, 3], order = [2, 1])
160 matC = triC
161 call disp%skip()
162 call disp%show("matA")
163 call disp%show( matA , format = cform )
164 call disp%show("matC")
165 call disp%show( matC , format = cform )
166 call disp%show("alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roffA = 0; coffA = 0; roffC = 0; coffC = 0;")
167 alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roffA = 0; coffA = 0; roffC = 0; coffC = 0;
168 call disp%show("call setMatUpdate(matC, symmetric, uppDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?SYRK contiguous interface.")
169 call setMatUpdate(matC, symmetric, uppDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?SYRK contiguous interface.
170 call disp%show("matC")
171 call disp%show( matC , format = cform )
172 call disp%show("matC - RefC")
173 call disp%show( matC - RefC , format = cform )
174 call disp%skip()
175 matC = triC
176 call disp%show("call setMatUpdate(matC(1:3, 1:3), symmetric, uppDia, matA(1:3, 1:5), nothing, alpha, beta) ! BLAS 3 ?SYRK simplified interface.")
177 call setMatUpdate(matC(1:3, 1:3), symmetric, uppDia, matA(1:3, 1:5), nothing, alpha, beta) ! BLAS 3 ?SYRK simplified interface.
178 call disp%show("matC")
179 call disp%show( matC , format = cform )
180 call disp%show("matC - RefC")
181 call disp%show( matC - RefC , format = cform )
182 call disp%skip()
183
184 end block
185
186 call disp%skip
187 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
188 call disp%show("!Complex hermitian matrix update.")
189 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
190 call disp%skip
191
192 block
193
194 complex(CKG) :: alpha, beta
195 complex(CKG), allocatable :: RefC(:,:), triC(:,:), matC(:,:), matA(:,:), VecX(:), VecY(:)
196
197 matA = reshape( [ complex(CKG) :: (2.0, 0.0), (3.0, 2.0), (4.0, 1.0) &
198 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
199 , (1.0, 3.0), (2.0, 1.0), (6.0, 0.0) &
200 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
201 , (1.0, 9.0), (3.0, 0.0), (6.0, 7.0) &
202 ], shape = [5, 3], order = [2, 1])
203 triC = reshape( [ complex(CKG) :: (6.0, DUM), CMPLX_DUMM, CMPLX_DUMM &
204 , (3.0, 4.0), (10.0, DUM), CMPLX_DUMM &
205 , (9.0, 1.0), (12.0, 2.0), (3.0, DUM) &
206 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
207 ], shape = [4, 3], order = [2, 1])
208 RefC = reshape( [ complex(CKG) :: (138.0, 0.0), CMPLX_DUMM, CMPLX_DUMM &
209 , (65.0, 80.0), (165.0, 0.0), CMPLX_DUMM &
210 , (134.0, 46.0), (88.0, -88.0), (199.0, 0.0) &
211 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
212 ], shape = [4, 3], order = [2, 1])
213 matC = triC
214 call disp%skip()
215 call disp%show("matA")
216 call disp%show( matA , format = cform )
217 call disp%show("matC")
218 call disp%show( matC , format = cform )
219 call disp%show("alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roffC = 0; coffC = 0; roffA = 0; coffA = 0;")
220 alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roffC = 0; coffC = 0; roffA = 0; coffA = 0;
221 call disp%show("call setMatUpdate(matC, hermitian, lowDia, matA, trans, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?HERK contiguous interface.")
222 call setMatUpdate(matC, hermitian, lowDia, matA, trans, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?HERK contiguous interface.
223 call disp%show("matC")
224 call disp%show( matC , format = cform )
225 call disp%show("matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.")
226 call disp%show( matC - RefC , format = cform )
227 call disp%skip()
228 matC = triC
229 call disp%show("call setMatUpdate(matC(1:3, 1:3), hermitian, lowDia, matA(1:5, 1:3), trans) ! BLAS 3 ?HERK simplified interface.")
230 call setMatUpdate(matC(1:3, 1:3), hermitian, lowDia, matA(1:5, 1:3), trans) ! BLAS 3 ?HERK simplified interface.
231 call disp%show("matC")
232 call disp%show( matC , format = cform )
233 call disp%show("matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.")
234 call disp%show( matC - RefC , format = cform )
235 call disp%skip()
236
237
238 matA = reshape( [ complex(CKG) :: (2.0, 0.0), (3.0, 2.0), (4.0, 1.0) &
239 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
240 , (1.0, 3.0), (2.0, 1.0), (6.0, 0.0) &
241 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
242 , (1.0, 9.0), (3.0, 0.0), (6.0, 7.0) &
243 ], shape = [3, 5], order = [1, 2])
244 matA = conjg(matA)
245 triC = reshape( [ complex(CKG) :: (6.0, DUM), CMPLX_DUMM, CMPLX_DUMM &
246 , (3.0, 4.0), (10.0, DUM), CMPLX_DUMM &
247 , (9.0, 1.0), (12.0, 2.0), (3.0, DUM) &
248 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
249 ], shape = [4, 3], order = [2, 1])
250 RefC = reshape( [ complex(CKG) :: (138.0, 0.0), CMPLX_DUMM, CMPLX_DUMM &
251 , (65.0, 80.0), (165.0, 0.0), CMPLX_DUMM &
252 , (134.0, 46.0), (88.0, -88.0), (199.0, 0.0) &
253 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
254 ], shape = [4, 3], order = [2, 1])
255 matC = triC
256 call disp%skip()
257 call disp%show("matA")
258 call disp%show( matA , format = cform )
259 call disp%show("matC")
260 call disp%show( matC , format = cform )
261 call disp%show("alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roffC = 0; coffC = 0; roffA = 0; coffA = 0;")
262 alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roffC = 0; coffC = 0; roffA = 0; coffA = 0;
263 call disp%show("call setMatUpdate(matC, hermitian, lowDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?HERK contiguous interface.")
264 call setMatUpdate(matC, hermitian, lowDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?HERK contiguous interface.
265 call disp%show("matC")
266 call disp%show( matC , format = cform )
267 call disp%show("matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.")
268 call disp%show( matC - RefC , format = cform )
269 call disp%skip()
270 matC = triC
271 call disp%show("call setMatUpdate(matC(1:3, 1:3), hermitian, lowDia, matA(1:3, 1:5), nothing) ! BLAS 3 ?HERK simplified interface.")
272 call setMatUpdate(matC(1:3, 1:3), hermitian, lowDia, matA(1:3, 1:5), nothing) ! BLAS 3 ?HERK simplified interface.
273 call disp%show("matC")
274 call disp%show( matC , format = cform )
275 call disp%show("matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.")
276 call disp%show( matC - RefC , format = cform )
277 call disp%skip()
278
279
280 matA = reshape( [ complex(CKG) :: (2.0, 0.0), (3.0, 2.0), (4.0, 1.0) &
281 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
282 , (1.0, 3.0), (2.0, 1.0), (6.0, 0.0) &
283 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
284 , (1.0, 9.0), (3.0, 0.0), (6.0, 7.0) &
285 ], shape = [3, 5])
286 matA = conjg(matA)
287 triC = reshape( [ complex(CKG) :: (6.0, DUM), CMPLX_DUMM, CMPLX_DUMM &
288 , (3.0, 4.0), (10.0, DUM), CMPLX_DUMM &
289 , (9.0, 1.0), (12.0, 2.0), (3.0, DUM) &
290 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
291 ], shape = [3, 4])
292 RefC = reshape( [ complex(CKG) :: (138., .0), (65., -72.), (134., -44.), CMPLX_DUMM &
293 , CMPLX_DUMM, (165.,0.), (88., 92.), CMPLX_DUMM &
294 , CMPLX_DUMM, CMPLX_DUMM, (199.,0.), CMPLX_DUMM &
295 ], shape = [3, 4], order = [2, 1])
296 matC = triC
297 call disp%skip()
298 call disp%show("matA")
299 call disp%show( matA , format = cform )
300 call disp%show("matC")
301 call disp%show( matC , format = cform )
302 call disp%show("alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roffA = 0; coffA = 0; roffC = 0; coffC = 0;")
303 alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roffA = 0; coffA = 0; roffC = 0; coffC = 0;
304 call disp%show("call setMatUpdate(matC, hermitian, uppDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?HERK contiguous interface.")
305 call setMatUpdate(matC, hermitian, uppDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?HERK contiguous interface.
306 call disp%show("matC")
307 call disp%show( matC , format = cform )
308 call disp%show("matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.")
309 call disp%show( matC - RefC , format = cform )
310 call disp%skip()
311 matC = triC
312 call disp%show("call setMatUpdate(matC(1:3, 1:3), hermitian, uppDia, matA(1:3, 1:5), nothing, alpha, beta) ! BLAS 3 ?HERK simplified interface.")
313 call setMatUpdate(matC(1:3, 1:3), hermitian, uppDia, matA(1:3, 1:5), nothing, alpha, beta) ! BLAS 3 ?HERK simplified interface.
314 call disp%show("matC")
315 call disp%show( matC , format = cform )
316 call disp%show("matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.")
317 call disp%show( matC - RefC , format = cform )
318 call disp%skip()
319
320 end block
321
322end program example
Generate and return a generic or type/kind-specific IO format with the requested specifications that ...
Definition: pm_io.F90:18485
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11726
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11508
This module contains classes and procedures for input/output (IO) or generic display operations on st...
Definition: pm_io.F90:252
type(display_type) disp
This is a scalar module variable an object of type display_type for general display.
Definition: pm_io.F90:11393
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
integer, parameter LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
Definition: pm_kind.F90:541
integer, parameter CKS
The single-precision complex kind in Fortran mode. On most platforms, this is a 32-bit real kind.
Definition: pm_kind.F90:570
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
Definition: pm_kind.F90:540
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
Definition: pm_kind.F90:539
integer, parameter RKS
The single-precision real kind in Fortran mode. On most platforms, this is an 32-bit real kind.
Definition: pm_kind.F90:567
Generate and return an object of type display_type.
Definition: pm_io.F90:10282

Example Unix compile command via Intel ifort compiler
1#!/usr/bin/env sh
2rm main.exe
3ifort -fpp -standard-semantics -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example Windows Batch compile command via Intel ifort compiler
1del main.exe
2set PATH=..\..\..\lib;%PATH%
3ifort /fpp /standard-semantics /O3 /I:..\..\..\include main.F90 ..\..\..\lib\libparamonte*.lib /exe:main.exe
4main.exe

Example Unix / MinGW compile command via GNU gfortran compiler
1#!/usr/bin/env sh
2rm main.exe
3gfortran -cpp -ffree-line-length-none -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example output
1
2matA
3+0.00000000, +8.00000000
4+1.00000000, +9.00000000
5+2.00000000, +10.0000000
6+3.00000000, +11.0000000
7+4.00000000, +12.0000000
8+5.00000000, +13.0000000
9+6.00000000, +14.0000000
10+7.00000000, +15.0000000
11-0.340282347E+39, -0.340282347E+39
12matC
13+0.00000000, +1.00000000, +3.00000000, +6.00000000, +10.0000000, +15.0000000, +21.0000000, +28.0000000
14-0.340282347E+39, +2.00000000, +4.00000000, +7.00000000, +11.0000000, +16.0000000, +22.0000000, +29.0000000
15-0.340282347E+39, -0.340282347E+39, +5.00000000, +8.00000000, +12.0000000, +17.0000000, +23.0000000, +30.0000000
16-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +9.00000000, +13.0000000, +18.0000000, +24.0000000, +31.0000000
17-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +14.0000000, +19.0000000, +25.0000000, +32.0000000
18-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +20.0000000, +26.0000000, +33.0000000
19-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +27.0000000, +34.0000000
20-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +35.0000000
21-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
22-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
23alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 2; roffC = 0; coffC = 0; roffA = 0; coffA = 0;
24call setMatUpdate(matC, symmetric, uppDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?SYRK contiguous interface.
25matC
26+64.0000000, +73.0000000, +83.0000000, +94.0000000, +106.000000, +119.000000, +133.000000, +148.000000
27-0.340282347E+39, +84.0000000, +96.0000000, +109.000000, +123.000000, +138.000000, +154.000000, +171.000000
28-0.340282347E+39, -0.340282347E+39, +109.000000, +124.000000, +140.000000, +157.000000, +175.000000, +194.000000
29-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +139.000000, +157.000000, +176.000000, +196.000000, +217.000000
30-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +174.000000, +195.000000, +217.000000, +240.000000
31-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +214.000000, +238.000000, +263.000000
32-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +259.000000, +286.000000
33-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +309.000000
34-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
35-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
36matC - RefC
37+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
38+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
39+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
40+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
41+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
42+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
43+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
44+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
45+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
46+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
47
48call setMatUpdate(matC(1:8,1:8), symmetric, uppDia, matA(1:8,1:2), nothing) ! BLAS 3 ?SYRK simplified interface.
49matC
50+64.0000000, +73.0000000, +83.0000000, +94.0000000, +106.000000, +119.000000, +133.000000, +148.000000
51-0.340282347E+39, +84.0000000, +96.0000000, +109.000000, +123.000000, +138.000000, +154.000000, +171.000000
52-0.340282347E+39, -0.340282347E+39, +109.000000, +124.000000, +140.000000, +157.000000, +175.000000, +194.000000
53-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +139.000000, +157.000000, +176.000000, +196.000000, +217.000000
54-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +174.000000, +195.000000, +217.000000, +240.000000
55-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +214.000000, +238.000000, +263.000000
56-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +259.000000, +286.000000
57-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, +309.000000
58-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
59-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
60matC - RefC
61+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
62+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
63+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
64+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
65+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
66+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
67+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
68+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
69+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
70+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
71
72
73matA
74+0.00000000, +3.00000000, +6.00000000, +9.00000000, +12.0000000, +15.0000000, +18.0000000, +21.0000000
75+1.00000000, +4.00000000, +7.00000000, +10.0000000, +13.0000000, +16.0000000, +19.0000000, +22.0000000
76+2.00000000, +5.00000000, +8.00000000, +11.0000000, +14.0000000, +17.0000000, +20.0000000, +23.0000000
77-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
78matC
79+0.00000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
80+1.00000000, +8.00000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
81+2.00000000, +9.00000000, +15.0000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
82+3.00000000, +10.0000000, +16.0000000, +21.0000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
83+4.00000000, +11.0000000, +17.0000000, +22.0000000, +26.0000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
84+5.00000000, +12.0000000, +18.0000000, +23.0000000, +27.0000000, +30.0000000, -0.340282347E+39, -0.340282347E+39
85+6.00000000, +13.0000000, +19.0000000, +24.0000000, +28.0000000, +31.0000000, +33.0000000, -0.340282347E+39
86+7.00000000, +14.0000000, +20.0000000, +25.0000000, +29.0000000, +32.0000000, +34.0000000, +35.0000000
87alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 3; roffA = 0; coffA = 0; roffC = 0; coffC = 0;
88call setMatUpdate(matC, symmetric, lowDia, matA, trans, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?SYRK contiguous interface.
89matC
90+5.00000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
91+15.0000000, +58.0000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
92+25.0000000, +95.0000000, +164.000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
93+35.0000000, +132.000000, +228.000000, +323.000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
94+45.0000000, +169.000000, +292.000000, +414.000000, +535.000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
95+55.0000000, +206.000000, +356.000000, +505.000000, +653.000000, +800.000000, -0.340282347E+39, -0.340282347E+39
96+65.0000000, +243.000000, +420.000000, +596.000000, +771.000000, +945.000000, +1118.00000, -0.340282347E+39
97+75.0000000, +280.000000, +484.000000, +687.000000, +889.000000, +1090.00000, +1290.00000, +1489.00000
98matC - RefC
99+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
100+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
101+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
102+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
103+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
104+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
105+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
106+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
107
108call setMatUpdate(matC(1:8, 1:8), symmetric, lowDia, matA(1:3, 1:8), trans) ! BLAS 3 ?SYRK simplified interface.
109matC
110+5.00000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
111+15.0000000, +58.0000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
112+25.0000000, +95.0000000, +164.000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
113+35.0000000, +132.000000, +228.000000, +323.000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
114+45.0000000, +169.000000, +292.000000, +414.000000, +535.000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
115+55.0000000, +206.000000, +356.000000, +505.000000, +653.000000, +800.000000, -0.340282347E+39, -0.340282347E+39
116+65.0000000, +243.000000, +420.000000, +596.000000, +771.000000, +945.000000, +1118.00000, -0.340282347E+39
117+75.0000000, +280.000000, +484.000000, +687.000000, +889.000000, +1090.00000, +1290.00000, +1489.00000
118matC - RefC
119+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
120+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
121+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
122+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
123+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
124+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
125+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
126+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
127
128
129!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
130!Complex Symmetric matrix update.
131!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
132
133
134matA
135( +2.000000, +0.000000), ( +3.000000, +2.000000), ( +4.000000, +1.000000), ( +1.000000, +7.000000), ( +0.000000, +0.000000)
136( +3.000000, +3.000000), ( +8.000000, +0.000000), ( +2.000000, +5.000000), ( +2.000000, +4.000000), ( +1.000000, +2.000000)
137( +1.000000, +3.000000), ( +2.000000, +1.000000), ( +6.000000, +0.000000), ( +3.000000, +2.000000), ( +2.000000, +2.000000)
138matC
139( +2.000000, +1.000000), ( +1.000000, +9.000000), ( +4.000000, +5.000000)
140(***************, ***************), ( +3.000000, +1.000000), ( +6.000000, +7.000000)
141(***************, ***************), (***************, ***************), ( +8.000000, +1.000000)
142(***************, ***************), (***************, ***************), (***************, ***************)
143alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roffA = 0; coffA = 0; roffC = 0; coffC = 0;
144call setMatUpdate(matC, symmetric, uppDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?SYRK contiguous interface.
145matC
146( -57.000000, +13.000000), ( -63.000000, +79.000000), ( -24.000000, +70.000000)
147(***************, ***************), ( -28.000000, +90.000000), ( -55.000000, +103.000000)
148(***************, ***************), (***************, ***************), ( +13.000000, +75.000000)
149(***************, ***************), (***************, ***************), (***************, ***************)
150matC - RefC
151( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
152( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
153( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
154( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
155
156call setMatUpdate(matC(1:3, 1:3), symmetric, uppDia, matA(1:3, 1:5), nothing, alpha, beta) ! BLAS 3 ?SYRK simplified interface.
157matC
158( -57.000000, +13.000000), ( -63.000000, +79.000000), ( -24.000000, +70.000000)
159(***************, ***************), ( -28.000000, +90.000000), ( -55.000000, +103.000000)
160(***************, ***************), (***************, ***************), ( +13.000000, +75.000000)
161(***************, ***************), (***************, ***************), (***************, ***************)
162matC - RefC
163( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
164( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
165( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
166( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
167
168
169!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
170!Complex hermitian matrix update.
171!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
172
173
174matA
175( +2.000000, +0.000000), ( +3.000000, +2.000000), ( +4.000000, +1.000000)
176( +3.000000, +3.000000), ( +8.000000, +0.000000), ( +2.000000, +5.000000)
177( +1.000000, +3.000000), ( +2.000000, +1.000000), ( +6.000000, +0.000000)
178( +3.000000, +3.000000), ( +8.000000, +0.000000), ( +2.000000, +5.000000)
179( +1.000000, +9.000000), ( +3.000000, +0.000000), ( +6.000000, +7.000000)
180matC
181( +6.000000, ***************), (***************, ***************), (***************, ***************)
182( +3.000000, +4.000000), ( +10.000000, ***************), (***************, ***************)
183( +9.000000, +1.000000), ( +12.000000, +2.000000), ( +3.000000, ***************)
184(***************, ***************), (***************, ***************), (***************, ***************)
185alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roffC = 0; coffC = 0; roffA = 0; coffA = 0;
186call setMatUpdate(matC, hermitian, lowDia, matA, trans, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?HERK contiguous interface.
187matC
188( +138.000000, +0.000000), (***************, ***************), (***************, ***************)
189( +65.000000, +80.000000), ( +165.000000, +0.000000), (***************, ***************)
190( +134.000000, +46.000000), ( +88.000000, -88.000000), ( +199.000000, +0.000000)
191(***************, ***************), (***************, ***************), (***************, ***************)
192matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.
193( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
194( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
195( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
196( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
197
198call setMatUpdate(matC(1:3, 1:3), hermitian, lowDia, matA(1:5, 1:3), trans) ! BLAS 3 ?HERK simplified interface.
199matC
200( +138.000000, +0.000000), (***************, ***************), (***************, ***************)
201( +65.000000, +80.000000), ( +165.000000, +0.000000), (***************, ***************)
202( +134.000000, +46.000000), ( +88.000000, -88.000000), ( +199.000000, +0.000000)
203(***************, ***************), (***************, ***************), (***************, ***************)
204matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.
205( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
206( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
207( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
208( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
209
210
211matA
212( +2.000000, -0.000000), ( +3.000000, -3.000000), ( +1.000000, -3.000000), ( +3.000000, -3.000000), ( +1.000000, -9.000000)
213( +3.000000, -2.000000), ( +8.000000, -0.000000), ( +2.000000, -1.000000), ( +8.000000, -0.000000), ( +3.000000, -0.000000)
214( +4.000000, -1.000000), ( +2.000000, -5.000000), ( +6.000000, -0.000000), ( +2.000000, -5.000000), ( +6.000000, -7.000000)
215matC
216( +6.000000, ***************), (***************, ***************), (***************, ***************)
217( +3.000000, +4.000000), ( +10.000000, ***************), (***************, ***************)
218( +9.000000, +1.000000), ( +12.000000, +2.000000), ( +3.000000, ***************)
219(***************, ***************), (***************, ***************), (***************, ***************)
220alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roffC = 0; coffC = 0; roffA = 0; coffA = 0;
221call setMatUpdate(matC, hermitian, lowDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?HERK contiguous interface.
222matC
223( +138.000000, +0.000000), (***************, ***************), (***************, ***************)
224( +65.000000, +80.000000), ( +165.000000, +0.000000), (***************, ***************)
225( +134.000000, +46.000000), ( +88.000000, -88.000000), ( +199.000000, +0.000000)
226(***************, ***************), (***************, ***************), (***************, ***************)
227matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.
228( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
229( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
230( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
231( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
232
233call setMatUpdate(matC(1:3, 1:3), hermitian, lowDia, matA(1:3, 1:5), nothing) ! BLAS 3 ?HERK simplified interface.
234matC
235( +138.000000, +0.000000), (***************, ***************), (***************, ***************)
236( +65.000000, +80.000000), ( +165.000000, +0.000000), (***************, ***************)
237( +134.000000, +46.000000), ( +88.000000, -88.000000), ( +199.000000, +0.000000)
238(***************, ***************), (***************, ***************), (***************, ***************)
239matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.
240( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
241( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
242( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
243( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
244
245
246matA
247( +2.000000, -0.000000), ( +3.000000, -3.000000), ( +1.000000, -3.000000), ( +3.000000, -3.000000), ( +1.000000, -9.000000)
248( +3.000000, -2.000000), ( +8.000000, -0.000000), ( +2.000000, -1.000000), ( +8.000000, -0.000000), ( +3.000000, -0.000000)
249( +4.000000, -1.000000), ( +2.000000, -5.000000), ( +6.000000, -0.000000), ( +2.000000, -5.000000), ( +6.000000, -7.000000)
250matC
251( +6.000000, ***************), ( +3.000000, +4.000000), ( +9.000000, +1.000000), (***************, ***************)
252(***************, ***************), ( +10.000000, ***************), ( +12.000000, +2.000000), (***************, ***************)
253(***************, ***************), (***************, ***************), ( +3.000000, ***************), (***************, ***************)
254alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roffA = 0; coffA = 0; roffC = 0; coffC = 0;
255call setMatUpdate(matC, hermitian, uppDia, matA, nothing, alpha, beta, ndim, ndum, roffC, coffC, roffA, coffA) ! BLAS 3 ?HERK contiguous interface.
256matC
257( +138.000000, +0.000000), ( +65.000000, -72.000000), ( +134.000000, -44.000000), (***************, ***************)
258(***************, ***************), ( +165.000000, +0.000000), ( +88.000000, +92.000000), (***************, ***************)
259(***************, ***************), (***************, ***************), ( +199.000000, +0.000000), (***************, ***************)
260matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.
261( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
262( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
263( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
264
265call setMatUpdate(matC(1:3, 1:3), hermitian, uppDia, matA(1:3, 1:5), nothing, alpha, beta) ! BLAS 3 ?HERK simplified interface.
266matC
267( +138.000000, +0.000000), ( +65.000000, -72.000000), ( +134.000000, -44.000000), (***************, ***************)
268(***************, ***************), ( +165.000000, +0.000000), ( +88.000000, +92.000000), (***************, ***************)
269(***************, ***************), (***************, ***************), ( +199.000000, +0.000000), (***************, ***************)
270matC - RefC ! Note the imaginary parts of the diagonals are set to zero on output.
271( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
272( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
273( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
274
275
Test:
test_pm_matrixUpdate
Todo:
Normal Priority: The input shape and offset arguments can be made optional by adding new interfaces to the module.
This should be done only after performing relevant benchmarks with the current interface to gauge whether the extension of this module is worth the effort.


Final Remarks


If you believe this algorithm or its documentation can be improved, we appreciate your contribution and help to edit this page's documentation and source file on GitHub.
For details on the naming abbreviations, see this page.
For details on the naming conventions, see this page.
This software is distributed under the MIT license with additional terms outlined below.

  1. If you use any parts or concepts from this library to any extent, please acknowledge the usage by citing the relevant publications of the ParaMonte library.
  2. If you regenerate any parts/ideas from this library in a programming environment other than those currently supported by this ParaMonte library (i.e., other than C, C++, Fortran, MATLAB, Python, R), please also ask the end users to cite this original ParaMonte library.

This software is available to the public under a highly permissive license.
Help us justify its continued development and maintenance by acknowledging its benefit to society, distributing it, and contributing to it.

Author:
Fatemeh Bagheri, Tuesday 08:49 PM, August 10, 2021, Dallas, TX Amir Shahmoradi, Friday 1:54 AM, April 21, 2017, Institute for Computational Engineering and Sciences (ICES), The University of Texas, Austin, TX

Definition at line 302 of file pm_matrixUpdate.F90.


The documentation for this interface was generated from the following file: