ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_matrixUpdate::setMatUpdateTriang 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 \(\texttt{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 \(\texttt{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*} & \texttt{mat} \leftarrow \alpha \texttt{matA}~~ \texttt{matA}^T + \beta \texttt{mat} &&\text{if operationA is missing.} \\ & \texttt{mat} \leftarrow \alpha \texttt{matA}^T \texttt{matA}~~ + \beta \texttt{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*} & \texttt{mat} \leftarrow \alpha \texttt{matA}~~ \texttt{matA}^H + \beta \texttt{mat} &&\text{if operationA is missing.} \\ & \texttt{mat} \leftarrow \alpha \texttt{matA}^H \texttt{matA}~~ + \beta \texttt{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*} & \texttt{mat} \leftarrow \beta \texttt{matA}~~ \texttt{matA}^T + \beta \texttt{mat} &&\text{if operationA is missing.} \\ & \texttt{mat} \leftarrow \beta \texttt{matA}^T \texttt{matA}~~ + \beta \texttt{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*} & \texttt{mat} \leftarrow \beta \texttt{matA}~~ \texttt{matA}^H + \beta \texttt{mat} &&\text{if operationA is missing.} \\ & \texttt{mat} \leftarrow \beta \texttt{matA}^H \texttt{matA}~~ + \beta \texttt{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]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.
[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.


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 setMatUpdateTriang(matA, mat, subset, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
call setMatUpdateTriang(matA, operationA, mat, subset, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
call setMatUpdateTriang(matA, operationA, 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 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 setMatUpdateTriang 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: transSymm, transHerm
9 use pm_io, only: getFormat
10 use pm_io, only: display_type
11
12 implicit none
13
14 type(display_type) :: disp
15
16 real(RKG) , parameter :: DUM = -huge(0._RKG)
17 complex(CKG) , parameter :: CMPLX_DUMM = cmplx(-huge(0._CKG), -huge(0._CKG), CKG)
18 integer(IK) :: ndim, ndum, roff, coff, roffA, coffA
19 character(:, SK), allocatable :: cform
20 cform = getFormat([cmplx(0., 0., CKG)], ed = SK_'f', signed = .true.)
21
22 disp = display_type(file = "main.out.F90")
23
24 block
25
26 real(RKG) :: alpha, beta
27 real(RKG), allocatable :: ref(:,:), mat(:,:), matA(:,:), VecX(:), VecY(:)
28
29 matA = reshape( [ real(RKG) :: 0.0, 8.0 &
30 , 1.0, 9.0 &
31 , 2.0, 10.0 &
32 , 3.0, 11.0 &
33 , 4.0, 12.0 &
34 , 5.0, 13.0 &
35 , 6.0, 14.0 &
36 , 7.0, 15.0 &
37 , DUM, DUM &
38 ], shape = [9, 2], order = [2, 1])
39 mat = reshape( [ real(RKG) :: 0.0, 1.0, 3.0, 6.0, 10.0, 15.0, 21.0, 28.0 &
40 , DUM, 2.0, 4.0, 7.0, 11.0, 16.0, 22.0, 29.0 &
41 , DUM, DUM, 5.0, 8.0, 12.0, 17.0, 23.0, 30.0 &
42 , DUM, DUM, DUM, 9.0, 13.0, 18.0, 24.0, 31.0 &
43 , DUM, DUM, DUM, DUM, 14.0, 19.0, 25.0, 32.0 &
44 , DUM, DUM, DUM, DUM, DUM, 20.0, 26.0, 33.0 &
45 , DUM, DUM, DUM, DUM, DUM, DUM, 27.0, 34.0 &
46 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, 35.0 &
47 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
48 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
49 ], shape = [10, 8], order = [2, 1])
50 ref = reshape( [ real(RKG) :: 64.0, 73.0, 83.0, 94.0, 106.0, 119.0, 133.0, 148.0 &
51 , DUM, 84.0, 96.0, 109.0, 123.0, 138.0, 154.0, 171.0 &
52 , DUM, DUM, 109.0, 124.0, 140.0, 157.0, 175.0, 194.0 &
53 , DUM, DUM, DUM, 139.0, 157.0, 176.0, 196.0, 217.0 &
54 , DUM, DUM, DUM, DUM, 174.0, 195.0, 217.0, 240.0 &
55 , DUM, DUM, DUM, DUM, DUM, 214.0, 238.0, 263.0 &
56 , DUM, DUM, DUM, DUM, DUM, DUM, 259.0, 286.0 &
57 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, 309.0 &
58 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
59 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
60 ], shape = [10, 8], order = [2, 1])
61 call disp%skip()
62 call disp%show("matA")
63 call disp%show( matA )
64 call disp%show("mat")
65 call disp%show( mat )
66 call disp%show("alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 2; roff = 0; coff = 0; roffA = 0; coffA = 0;")
67 alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 2; roff = 0; coff = 0; roffA = 0; coffA = 0;
68 call disp%show("call setMatUpdateTriang(mat, uppDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)")
69 call setMatUpdateTriang(mat, uppDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
70 call disp%show("mat")
71 call disp%show( mat )
72 call disp%show("mat - ref")
73 call disp%show( mat - ref )
74 call disp%skip()
75
76 matA = reshape( [ real(RKG) :: 0.0, 3.0, 6.0, 9.0, 12.0, 15.0, 18.0, 21.0 &
77 , 1.0, 4.0, 7.0, 10.0, 13.0, 16.0, 19.0, 22.0 &
78 , 2.0, 5.0, 8.0, 11.0, 14.0, 17.0, 20.0, 23.0 &
79 , DUM, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
80 ], shape = [4, 8], order = [2, 1])
81 mat = reshape( [ real(RKG) :: 0.0, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
82 , 1.0, 8.0, DUM, DUM, DUM, DUM, DUM, DUM &
83 , 2.0, 9.0, 15.0, DUM, DUM, DUM, DUM, DUM &
84 , 3.0, 10.0, 16.0, 21.0, DUM, DUM, DUM, DUM &
85 , 4.0, 11.0, 17.0, 22.0, 26.0, DUM, DUM, DUM &
86 , 5.0, 12.0, 18.0, 23.0, 27.0, 30.0, DUM, DUM &
87 , 6.0, 13.0, 19.0, 24.0, 28.0, 31.0, 33.0, DUM &
88 , 7.0, 14.0, 20.0, 25.0, 29.0, 32.0, 34.0, 35.0 &
89 ], shape = [8, 8], order = [2, 1])
90 ref = reshape( [ real(RKG) :: 5.0, DUM, DUM, DUM, DUM, DUM, DUM, DUM &
91 , 15.0, 58.0, DUM, DUM, DUM, DUM, DUM, DUM &
92 , 25.0, 95.0, 164.0, DUM, DUM, DUM, DUM, DUM &
93 , 35.0, 132.0, 228.0, 323.0, DUM, DUM, DUM, DUM &
94 , 45.0, 169.0, 292.0, 414.0, 535.0, DUM, DUM, DUM &
95 , 55.0, 206.0, 356.0, 505.0, 653.0, 800.0, DUM, DUM &
96 , 65.0, 243.0, 420.0, 596.0, 771.0, 945.0, 1118.0, DUM &
97 , 75.0, 280.0, 484.0, 687.0, 889.0, 1090.0, 1290.0, 1489.0 &
98 ], shape = [8, 8], order = [2, 1])
99 call disp%skip()
100 call disp%show("matA")
101 call disp%show( matA )
102 call disp%show("mat")
103 call disp%show( mat )
104 call disp%show("alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 3; roff = 0; coff = 0; roffA = 0; coffA = 0;")
105 alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 3; roff = 0; coff = 0; roffA = 0; coffA = 0;
106 call disp%show("call setMatUpdateTriang(mat, lowDia, matA, transSymm, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)")
107 call setMatUpdateTriang(mat, lowDia, matA, transSymm, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
108 call disp%show("mat")
109 call disp%show( mat )
110 call disp%show("mat - ref")
111 call disp%show( mat - ref )
112 call disp%skip()
113
114 end block
115
116 call disp%skip
117 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
118 call disp%show("!Complex Symmetric matrix update.")
119 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
120 call disp%skip
121
122 block
123
124 complex(CKG) :: alpha, beta
125 complex(CKG), allocatable :: ref(:,:), mat(:,:), matA(:,:), VecX(:), VecY(:)
126
127 matA = reshape( [ complex(CKG) :: (2.0, 0.0), (3.0, 2.0), (4.0, 1.0), (1.0, 7.0), (0.0, 0.0) &
128 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0), (2.0, 4.0), (1.0, 2.0) &
129 , (1.0, 3.0), (2.0, 1.0), (6.0, 0.0), (3.0, 2.0), (2.0, 2.0) &
130 ], shape = [3, 5], order = [2, 1])
131 mat = reshape( [ complex(CKG) :: (2.0, 1.0), (1.0, 9.0), (4.0, 5.0) &
132 , CMPLX_DUMM, (3.0, 1.0), (6.0, 7.0) &
133 , CMPLX_DUMM, CMPLX_DUMM, (8.0, 1.0) &
134 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
135 ], shape = [4, 3], order = [2, 1])
136 ref = reshape( [ complex(CKG) :: (-57.0, 13.0), (-63.0, 79.0), (-24.0, 70.0) &
137 , CMPLX_DUMM, (-28.0, 90.0), (-55.0, 103.0) &
138 , CMPLX_DUMM, CMPLX_DUMM, (13.0, 75.0) &
139 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
140 ], shape = [4, 3], order = [2, 1])
141 call disp%skip()
142 call disp%show("matA")
143 call disp%show( matA , format = cform )
144 call disp%show("mat")
145 call disp%show( mat , format = cform )
146 call disp%show("alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;")
147 alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;
148 call disp%show("call setMatUpdateTriang(mat, uppDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)")
149 call setMatUpdateTriang(mat, uppDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
150 call disp%show("mat")
151 call disp%show( mat , format = cform )
152 call disp%show("mat - ref")
153 call disp%show( mat - ref , format = cform )
154 call disp%skip()
155
156 end block
157
158 call disp%skip
159 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
160 call disp%show("!Complex hermitian matrix update.")
161 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
162 call disp%skip
163
164 block
165
166 real(CKG) :: alpha, beta
167 complex(CKG), allocatable :: ref(:,:), mat(:,:), matA(:,:), VecX(:), VecY(:)
168
169 matA = reshape( [ complex(CKG) :: (2.0, 0.0), (3.0, 2.0), (4.0, 1.0) &
170 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
171 , (1.0, 3.0), (2.0, 1.0), (6.0, 0.0) &
172 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
173 , (1.0, 9.0), (3.0, 0.0), (6.0, 7.0) &
174 ], shape = [5, 3], order = [2, 1])
175 mat = reshape( [ complex(CKG) :: (6.0, DUM), CMPLX_DUMM, CMPLX_DUMM &
176 , (3.0, 4.0), (10.0, DUM), CMPLX_DUMM &
177 , (9.0, 1.0), (12.0, 2.0), (3.0, DUM) &
178 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
179 ], shape = [4, 3], order = [2, 1])
180 ref = reshape( [ complex(CKG) :: (138.0, 0.0), CMPLX_DUMM, CMPLX_DUMM &
181 , (65.0, 80.0), (165.0, 0.0), CMPLX_DUMM &
182 , (134.0, 46.0), (88.0, -88.0), (199.0, 0.0) &
183 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
184 ], shape = [4, 3], order = [2, 1])
185 call disp%skip()
186 call disp%show("matA")
187 call disp%show( matA , format = cform )
188 call disp%show("mat")
189 call disp%show( mat , format = cform )
190 call disp%show("alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;")
191 alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;
192 call disp%show("call setMatUpdateTriang(mat, lowDia, matA, transHerm, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)")
193 call setMatUpdateTriang(mat, lowDia, matA, transHerm, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
194 call disp%show("mat")
195 call disp%show( mat , format = cform )
196 call disp%show("mat - ref ! Note the imaginary parts of the diagonals are set to zero on output.")
197 call disp%show( mat - ref , format = cform )
198 call disp%skip()
199
200
201 matA = reshape( [ complex(CKG) :: (2.0, 0.0), (3.0, 2.0), (4.0, 1.0) &
202 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
203 , (1.0, 3.0), (2.0, 1.0), (6.0, 0.0) &
204 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
205 , (1.0, 9.0), (3.0, 0.0), (6.0, 7.0) &
206 ], shape = [3, 5], order = [1, 2])
207 matA = conjg(matA)
208 mat = reshape( [ complex(CKG) :: (6.0, DUM), CMPLX_DUMM, CMPLX_DUMM &
209 , (3.0, 4.0), (10.0, DUM), CMPLX_DUMM &
210 , (9.0, 1.0), (12.0, 2.0), (3.0, DUM) &
211 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
212 ], shape = [4, 3], order = [2, 1])
213 ref = reshape( [ complex(CKG) :: (138.0, 0.0), CMPLX_DUMM, CMPLX_DUMM &
214 , (65.0, 80.0), (165.0, 0.0), CMPLX_DUMM &
215 , (134.0, 46.0), (88.0, -88.0), (199.0, 0.0) &
216 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
217 ], shape = [4, 3], order = [2, 1])
218 call disp%skip()
219 call disp%show("matA")
220 call disp%show( matA , format = cform )
221 call disp%show("mat")
222 call disp%show( mat , format = cform )
223 call disp%show("alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;")
224 alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;
225 call disp%show("call setMatUpdateTriang(mat, lowDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)")
226 call setMatUpdateTriang(mat, lowDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
227 call disp%show("mat")
228 call disp%show( mat , format = cform )
229 call disp%show("mat - ref ! Note the imaginary parts of the diagonals are set to zero on output.")
230 call disp%show( mat - ref , format = cform )
231 call disp%skip()
232
233
234 matA = reshape( [ complex(CKG) :: (2.0, 0.0), (3.0, 2.0), (4.0, 1.0) &
235 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
236 , (1.0, 3.0), (2.0, 1.0), (6.0, 0.0) &
237 , (3.0, 3.0), (8.0, 0.0), (2.0, 5.0) &
238 , (1.0, 9.0), (3.0, 0.0), (6.0, 7.0) &
239 ], shape = [3, 5])
240 matA = conjg(matA)
241 mat = reshape( [ complex(CKG) :: (6.0, DUM), CMPLX_DUMM, CMPLX_DUMM &
242 , (3.0, 4.0), (10.0, DUM), CMPLX_DUMM &
243 , (9.0, 1.0), (12.0, 2.0), (3.0, DUM) &
244 , CMPLX_DUMM, CMPLX_DUMM, CMPLX_DUMM &
245 ], shape = [3, 4])
246 ref = reshape( [ complex(CKG) :: (138., .0), (65., -72.), (134., -44.), CMPLX_DUMM &
247 , CMPLX_DUMM, (165.,0.), (88., 92.), CMPLX_DUMM &
248 , CMPLX_DUMM, CMPLX_DUMM, (199.,0.), CMPLX_DUMM &
249 ], shape = [3, 4], order = [2, 1])
250 call disp%skip()
251 call disp%show("matA")
252 call disp%show( matA , format = cform )
253 call disp%show("mat")
254 call disp%show( mat , format = cform )
255 call disp%show("alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;")
256 alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;
257 call disp%show("call setMatUpdateTriang(mat, uppDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)")
258 call setMatUpdateTriang(mat, uppDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
259 call disp%show("mat")
260 call disp%show( mat , format = cform )
261 call disp%show("mat - ref ! Note the imaginary parts of the diagonals are set to zero on output.")
262 call disp%show( mat - ref , format = cform )
263 call disp%skip()
264
265 end block
266
267end 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
12mat
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; roff = 0; coff = 0; roffA = 0; coffA = 0;
24call setMatUpdateTriang(mat, uppDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
25mat
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
36mat - ref
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
48
49matA
50+0.00000000, +3.00000000, +6.00000000, +9.00000000, +12.0000000, +15.0000000, +18.0000000, +21.0000000
51+1.00000000, +4.00000000, +7.00000000, +10.0000000, +13.0000000, +16.0000000, +19.0000000, +22.0000000
52+2.00000000, +5.00000000, +8.00000000, +11.0000000, +14.0000000, +17.0000000, +20.0000000, +23.0000000
53-0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
54mat
55+0.00000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
56+1.00000000, +8.00000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
57+2.00000000, +9.00000000, +15.0000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
58+3.00000000, +10.0000000, +16.0000000, +21.0000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
59+4.00000000, +11.0000000, +17.0000000, +22.0000000, +26.0000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
60+5.00000000, +12.0000000, +18.0000000, +23.0000000, +27.0000000, +30.0000000, -0.340282347E+39, -0.340282347E+39
61+6.00000000, +13.0000000, +19.0000000, +24.0000000, +28.0000000, +31.0000000, +33.0000000, -0.340282347E+39
62+7.00000000, +14.0000000, +20.0000000, +25.0000000, +29.0000000, +32.0000000, +34.0000000, +35.0000000
63alpha = 1._RKG; beta = 1._RKG; ndim = 8; ndum = 3; roff = 0; coff = 0; roffA = 0; coffA = 0;
64call setMatUpdateTriang(mat, lowDia, matA, transSymm, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
65mat
66+5.00000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
67+15.0000000, +58.0000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
68+25.0000000, +95.0000000, +164.000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
69+35.0000000, +132.000000, +228.000000, +323.000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
70+45.0000000, +169.000000, +292.000000, +414.000000, +535.000000, -0.340282347E+39, -0.340282347E+39, -0.340282347E+39
71+55.0000000, +206.000000, +356.000000, +505.000000, +653.000000, +800.000000, -0.340282347E+39, -0.340282347E+39
72+65.0000000, +243.000000, +420.000000, +596.000000, +771.000000, +945.000000, +1118.00000, -0.340282347E+39
73+75.0000000, +280.000000, +484.000000, +687.000000, +889.000000, +1090.00000, +1290.00000, +1489.00000
74mat - ref
75+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
76+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
77+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
78+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
79+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
80+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
81+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
82+0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000, +0.00000000
83
84
85!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
86!Complex Symmetric matrix update.
87!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
88
89
90matA
91( +2.000000, +0.000000), ( +3.000000, +2.000000), ( +4.000000, +1.000000), ( +1.000000, +7.000000), ( +0.000000, +0.000000)
92( +3.000000, +3.000000), ( +8.000000, +0.000000), ( +2.000000, +5.000000), ( +2.000000, +4.000000), ( +1.000000, +2.000000)
93( +1.000000, +3.000000), ( +2.000000, +1.000000), ( +6.000000, +0.000000), ( +3.000000, +2.000000), ( +2.000000, +2.000000)
94mat
95( +2.000000, +1.000000), ( +1.000000, +9.000000), ( +4.000000, +5.000000)
96(***************, ***************), ( +3.000000, +1.000000), ( +6.000000, +7.000000)
97(***************, ***************), (***************, ***************), ( +8.000000, +1.000000)
98(***************, ***************), (***************, ***************), (***************, ***************)
99alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;
100call setMatUpdateTriang(mat, uppDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
101mat
102( -57.000000, +13.000000), ( -63.000000, +79.000000), ( -24.000000, +70.000000)
103(***************, ***************), ( -28.000000, +90.000000), ( -55.000000, +103.000000)
104(***************, ***************), (***************, ***************), ( +13.000000, +75.000000)
105(***************, ***************), (***************, ***************), (***************, ***************)
106mat - ref
107( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
108( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
109( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
110( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
111
112
113!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
114!Complex hermitian matrix update.
115!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
116
117
118matA
119( +2.000000, +0.000000), ( +3.000000, +2.000000), ( +4.000000, +1.000000)
120( +3.000000, +3.000000), ( +8.000000, +0.000000), ( +2.000000, +5.000000)
121( +1.000000, +3.000000), ( +2.000000, +1.000000), ( +6.000000, +0.000000)
122( +3.000000, +3.000000), ( +8.000000, +0.000000), ( +2.000000, +5.000000)
123( +1.000000, +9.000000), ( +3.000000, +0.000000), ( +6.000000, +7.000000)
124mat
125( +6.000000, ***************), (***************, ***************), (***************, ***************)
126( +3.000000, +4.000000), ( +10.000000, ***************), (***************, ***************)
127( +9.000000, +1.000000), ( +12.000000, +2.000000), ( +3.000000, ***************)
128(***************, ***************), (***************, ***************), (***************, ***************)
129alpha = (1._CKG, 1._CKG); beta = (1._CKG, 1._CKG); ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;
130call setMatUpdateTriang(mat, lowDia, matA, transHerm, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
131mat
132( +138.000000, +0.000000), (***************, ***************), (***************, ***************)
133( +65.000000, +80.000000), ( +165.000000, +0.000000), (***************, ***************)
134( +134.000000, +46.000000), ( +88.000000, -88.000000), ( +199.000000, +0.000000)
135(***************, ***************), (***************, ***************), (***************, ***************)
136mat - ref ! Note the imaginary parts of the diagonals are set to zero on output.
137( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
138( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
139( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
140( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
141
142
143matA
144( +2.000000, -0.000000), ( +3.000000, -3.000000), ( +1.000000, -3.000000), ( +3.000000, -3.000000), ( +1.000000, -9.000000)
145( +3.000000, -2.000000), ( +8.000000, -0.000000), ( +2.000000, -1.000000), ( +8.000000, -0.000000), ( +3.000000, -0.000000)
146( +4.000000, -1.000000), ( +2.000000, -5.000000), ( +6.000000, -0.000000), ( +2.000000, -5.000000), ( +6.000000, -7.000000)
147mat
148( +6.000000, ***************), (***************, ***************), (***************, ***************)
149( +3.000000, +4.000000), ( +10.000000, ***************), (***************, ***************)
150( +9.000000, +1.000000), ( +12.000000, +2.000000), ( +3.000000, ***************)
151(***************, ***************), (***************, ***************), (***************, ***************)
152alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;
153call setMatUpdateTriang(mat, lowDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
154mat
155( +138.000000, +0.000000), (***************, ***************), (***************, ***************)
156( +65.000000, +80.000000), ( +165.000000, +0.000000), (***************, ***************)
157( +134.000000, +46.000000), ( +88.000000, -88.000000), ( +199.000000, +0.000000)
158(***************, ***************), (***************, ***************), (***************, ***************)
159mat - ref ! Note the imaginary parts of the diagonals are set to zero on output.
160( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
161( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
162( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
163( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
164
165
166matA
167( +2.000000, -0.000000), ( +3.000000, -3.000000), ( +1.000000, -3.000000), ( +3.000000, -3.000000), ( +1.000000, -9.000000)
168( +3.000000, -2.000000), ( +8.000000, -0.000000), ( +2.000000, -1.000000), ( +8.000000, -0.000000), ( +3.000000, -0.000000)
169( +4.000000, -1.000000), ( +2.000000, -5.000000), ( +6.000000, -0.000000), ( +2.000000, -5.000000), ( +6.000000, -7.000000)
170mat
171( +6.000000, ***************), ( +3.000000, +4.000000), ( +9.000000, +1.000000), (***************, ***************)
172(***************, ***************), ( +10.000000, ***************), ( +12.000000, +2.000000), (***************, ***************)
173(***************, ***************), (***************, ***************), ( +3.000000, ***************), (***************, ***************)
174alpha = 1._CKG; beta = 1._CKG; ndim = 3; ndum = 5; roff = 0; coff = 0; roffA = 0; coffA = 0;
175call setMatUpdateTriang(mat, uppDia, matA, alpha, beta, ndim, ndum, roff, coff, roffA, coffA)
176mat
177( +138.000000, +0.000000), ( +65.000000, -72.000000), ( +134.000000, -44.000000), (***************, ***************)
178(***************, ***************), ( +165.000000, +0.000000), ( +88.000000, +92.000000), (***************, ***************)
179(***************, ***************), (***************, ***************), ( +199.000000, +0.000000), (***************, ***************)
180mat - ref ! Note the imaginary parts of the diagonals are set to zero on output.
181( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
182( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
183( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000), ( +0.000000, +0.000000)
184
185
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 5095 of file pm_matrixUpdate.F90.


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