2
3
4
5
6
7
8
9 INTEGER IA, JA, M, N
10
11
12 INTEGER DESCA( * )
13 COMPLEX*16 A( * ), TAU( * ), WORK( * )
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122 INTEGER BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_,
123 $ LLD_, MB_, M_, NB_, N_, RSRC_
124 parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
125 $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
126 $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
127 COMPLEX*16 ZERO
128 parameter( zero = ( 0.0d+0, 0.0d+0 ) )
129
130
131 CHARACTER COLBTOP, ROWBTOP
132 INTEGER I, IACOL, IAROW, IB, ICOFF, ICTXT, IIA, IN,
133 $ IPT, IPV, IPW, JJA, JM1, JV, L, MYCOL, MYROW,
134 $ NPCOL, NPROW, NQ
135
136
137 INTEGER DESCV( DLEN_ )
138
139
143
144
145 INTEGER ICEIL, NUMROC
147
148
150
151
152
153
154
155 ictxt = desca( ctxt_ )
156 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
157
158
159
160 IF( n.LT.m )
161 $ RETURN
162
163 l = n - m
164 jm1 = ja +
min( m+1, n ) - 1
165 in =
min(
iceil( ia, desca( mb_ ) ) * desca( mb_ ), ia+m-1 )
166 icoff = mod( ja-1, desca( nb_ ) )
167 CALL infog2l( ia, ja, desca, nprow, npcol, myrow, mycol, iia, jja,
168 $ iarow, iacol )
169 nq =
numroc( n+icoff, desca( nb_ ), mycol, iacol, npcol )
170 ipv = 1
171 ipt = ipv + nq * desca( mb_ )
172 ipw = ipt + desca( mb_ ) * desca( mb_ )
173 CALL pb_topget( ictxt, 'Broadcast', 'Rowwise', rowbtop )
174 CALL pb_topget( ictxt, 'Broadcast', 'Columnwise', colbtop )
175 CALL pb_topset( ictxt, 'Broadcast', 'Rowwise', ' ' )
176 CALL pb_topset( ictxt, 'Broadcast', 'Columnwise', 'I-ring' )
177
178 CALL descset( descv, desca( mb_ ), n + icoff, desca( mb_ ),
179 $ desca( nb_ ), iarow, iacol, ictxt, desca( mb_ ) )
180
181
182
183 ib = in - ia + 1
184 jv = icoff + jm1 - ja + 1
185
186
187
188 CALL pzlarzt(
'Backward',
'Rowwise', l, ib, a, ia, jm1, desca,
189 $ tau, work( ipt ), work( ipw ) )
190
191
192
193 CALL pzlacpy(
'All', ib, l, a, ia, jm1, desca, work( ipv ), 1,
194 $ jv, descv )
195
196
197
198 CALL pzlacpy(
'Lower', ib-1, ib-1, a, ia+1, ja, desca,
199 $ work( ipv ), 1, icoff+1, descv )
200
201
202
203 CALL pzlaset(
'All', ib, l, zero, zero, a, ia, jm1, desca )
204 CALL pzlaset(
'Lower', ib-1, ib-1, zero, zero, a, ia+1, ja,
205 $ desca )
206
207
208
209 CALL pzlarzb(
'Right',
'Conjugate transpose',
'Backward',
210 $ 'Rowwise', in-ia+1, n, ib, l, work( ipv ), 1, jv,
211 $ descv, work( ipt ), a, ia, ja, desca, work( ipw ) )
212
213
214
215 CALL pzlacpy(
'Lower', ib-1, ib-1, work( ipv ), 1, icoff+1, descv,
216 $ a, ia+1, ja, desca )
217
218 descv( rsrc_ ) = mod( descv( rsrc_ ) + 1, nprow )
219
220
221
222 DO 10 i = in+1, ia+m-1, desca( mb_ )
223 ib =
min( ia+m-i, desca( mb_ ) )
224
225
226
227 CALL pzlarzt(
'Backward',
'Rowwise', l, ib, a, i, jm1, desca,
228 $ tau, work( ipt ), work( ipw ) )
229
230
231
232 CALL pzlacpy(
'All', ib, l, a, i, jm1, desca, work( ipv ), 1,
233 $ jv, descv )
234
235
236
237 CALL pzlacpy(
'Lower', ib-1, ib-1, a, i+1, ja+i-ia, desca,
238 $ work( ipv ), 1, icoff+1+i-ia, descv )
239
240
241
242 CALL pzlaset(
'All', ib, l, zero, zero, a, i, jm1, desca )
243 CALL pzlaset(
'Lower', ib-1, ib-1, zero, zero, a, i+1, ja+i-ia,
244 $ desca )
245
246
247
248 CALL pzlarzb(
'Right',
'Conjugate transpose',
'Backward',
249 $ 'Rowwise', i+ib-ia, n-i+ia, ib, l, work( ipv ),
250 $ 1, jv, descv, work( ipt ), a, ia, ja+i-ia, desca,
251 $ work( ipw ) )
252
253 CALL pzlacpy(
'Lower', ib-1, ib-1, work( ipv ), 1,
254 $ icoff+1+i-ia, descv, a, i+1, ja+i-ia, desca )
255
256 descv( rsrc_ ) = mod( descv( rsrc_ ) + 1, nprow )
257
258 10 CONTINUE
259
260 CALL pb_topset( ictxt, 'Broadcast', 'Rowwise', rowbtop )
261 CALL pb_topset( ictxt, 'Broadcast', 'Columnwise', colbtop )
262
263 RETURN
264
265
266
subroutine descset(desc, m, n, mb, nb, irsrc, icsrc, ictxt, lld)
integer function iceil(inum, idenom)
subroutine infog2l(grindx, gcindx, desc, nprow, npcol, myrow, mycol, lrindx, lcindx, rsrc, csrc)
integer function numroc(n, nb, iproc, isrcproc, nprocs)
subroutine pzlaset(uplo, m, n, alpha, beta, a, ia, ja, desca)
subroutine pzlacpy(uplo, m, n, a, ia, ja, desca, b, ib, jb, descb)
subroutine pzlarzb(side, trans, direct, storev, m, n, k, l, v, iv, jv, descv, t, c, ic, jc, descc, work)
subroutine pzlarzt(direct, storev, n, k, v, iv, jv, descv, tau, t, work)