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 ONE, ZERO
128 parameter( one = ( 1.0d+0, 0.0d+0 ),
129 $ zero = ( 0.0d+0, 0.0d+0 ) )
130
131
132 CHARACTER COLBTOP, ROWBTOP
133 INTEGER I, IACOL, IAROW, IB, ICOFF, ICTXT, IIA, IN,
134 $ IPT, IPV, IPW, JJA, JV, K, MYCOL, MYROW, NPCOL,
135 $ NPROW, NQ
136
137
138 INTEGER DESCV( DLEN_ )
139
140
144
145
146 INTEGER ICEIL, NUMROC
148
149
151
152
153
154
155
156 ictxt = desca( ctxt_ )
157 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
158
160 in =
min(
iceil( ia+m-k, desca( mb_ ) ) * desca( mb_ ), ia+m-1 )
161
162 icoff = mod( ja-1, desca( nb_ ) )
163 CALL infog2l( ia+m-k, ja, desca, nprow, npcol, myrow, mycol,
164 $ iia, jja, iarow, iacol )
165 nq =
numroc( n+icoff, desca( nb_ ), mycol, iacol, npcol )
166 ipv = 1
167 ipt = ipv + nq * desca( mb_ )
168 ipw = ipt + desca( mb_ ) * desca( mb_ )
169 CALL pb_topget( ictxt, 'Broadcast', 'Rowwise', rowbtop )
170 CALL pb_topget( ictxt, 'Broadcast', 'Columnwise', colbtop )
171 CALL pb_topset( ictxt, 'Broadcast', 'Rowwise', ' ' )
172 CALL pb_topset( ictxt, 'Broadcast', 'Columnwise', 'I-ring' )
173
174 CALL descset( descv, desca( mb_), n + icoff, desca( mb_ ),
175 $ desca( nb_ ), iarow, iacol, ictxt, desca( mb_ ) )
176
177
178
179 ib = in - ia - m + k + 1
180 jv = 1 + n - k + icoff
181
182
183
184 CALL pzlarft(
'Backward',
'Rowwise', n-m+in-ia+1, ib, a, ia+m-k,
185 $ ja, desca, tau, work( ipt ), work( ipw ) )
186
187
188
189 CALL pzlacpy(
'All', ib, n-m+in-ia+1, a, ia+m-k, ja, desca,
190 $ work( ipv ), 1, icoff+1, descv )
191 CALL pzlaset(
'Upper', ib, ib, zero, one, work( ipv ), 1, jv,
192 $ descv )
193
194
195
196
197 CALL pzlaset(
'All', ib, n-k, zero, zero, a, ia+m-k, ja,
198 $ desca )
199 CALL pzlaset(
'Lower', ib-1, ib, zero, zero, a, ia+m-k+1,
200 $ ja+n-k, desca )
201
202
203
204 CALL pzlarfb(
'Right',
'Conjugate transpose',
'Backward',
205 $ 'Rowwise', in-ia+1, n-m+in-ia+1, ib, work( ipv ), 1,
206 $ icoff+1, descv, work( ipt ), a, ia, ja, desca,
207 $ work( ipw ) )
208
209 descv( rsrc_ ) = mod( descv( rsrc_ ) + 1, nprow )
210
211
212
213 DO 10 i = in+1, ia+m-1, desca( mb_ )
214 ib =
min( ia+m-i, desca( mb_ ) )
215 jv = 1 + n - m + i - ia + icoff
216
217
218
219 CALL pzlarft(
'Backward',
'Rowwise', n-m+i+ib-ia, ib, a, i, ja,
220 $ desca, tau, work( ipt ), work( ipw ) )
221
222
223
224 CALL pzlacpy(
'All', ib, n-m+i+ib-ia, a, i, ja, desca,
225 $ work( ipv ), 1, icoff+1, descv )
226 CALL pzlaset(
'Upper', ib, ib, zero, one, work( ipv ), 1, jv,
227 $ descv )
228
229
230
231
232 CALL pzlaset(
'All', ib, n-m+i-ia, zero, zero, a, i, ja,
233 $ desca )
234 CALL pzlaset(
'Lower', ib-1, ib, zero, zero, a, i+1,
235 $ ja+n-m+i-ia, desca )
236
237
238
239 CALL pzlarfb(
'Right',
'Conjugate transpose',
'Backward',
240 $ 'Rowwise', i+ib-ia, n-m+i+ib-ia, ib, work( ipv ),
241 $ 1, icoff+1, descv, work( ipt ), a, ia, ja, desca,
242 $ work( ipw ) )
243
244 descv( rsrc_ ) = mod( descv( rsrc_ ) + 1, nprow )
245
246 10 CONTINUE
247
248 CALL pb_topset( ictxt, 'Broadcast', 'Rowwise', rowbtop )
249 CALL pb_topset( ictxt, 'Broadcast', 'Columnwise', colbtop )
250
251 RETURN
252
253
254
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 pzlarfb(side, trans, direct, storev, m, n, k, v, iv, jv, descv, t, c, ic, jc, descc, work)
subroutine pzlarft(direct, storev, n, k, v, iv, jv, descv, tau, t, work)