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
123
124 INTEGER BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_,
125 $ LLD_, MB_, M_, NB_, N_, RSRC_
126 parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
127 $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
128 $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
129 COMPLEX*16 ONE, ZERO
130 parameter( one = ( 1.0d+0, 0.0d+0 ),
131 $ zero = ( 0.0d+0, 0.0d+0 ) )
132
133
134 CHARACTER COLBTOP, ROWBTOP
135 INTEGER IACOL, IAROW, ICTXT, IIA, IPT, IPV, IPW, IROFF,
136 $ IV, J, JB, JJA, JN, K, MP, MYCOL, MYROW, NPCOL,
137 $ NPROW
138
139
140 INTEGER DESCV( DLEN_ )
141
142
146
147
148 INTEGER ICEIL, NUMROC
150
151
153
154
155
156
157
158 ictxt = desca( ctxt_ )
159 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
160
162 jn =
min(
iceil( ja+n-k, desca( nb_ ) ) * desca( nb_ ), ja+n-1 )
163
164 iroff = mod( ia-1, desca( mb_ ) )
165 CALL infog2l( ia, ja+n-k, desca, nprow, npcol, myrow, mycol,
166 $ iia, jja, iarow, iacol )
167 mp =
numroc( m+iroff, desca( mb_ ), myrow, iarow, nprow )
168 ipv = 1
169 ipt = ipv + mp * desca( nb_ )
170 ipw = ipt + desca( nb_ ) * desca( nb_ )
171 CALL pb_topget( ictxt, 'Broadcast', 'Rowwise', rowbtop )
172 CALL pb_topget( ictxt, 'Broadcast', 'Columnwise', colbtop )
173 CALL pb_topset( ictxt, 'Broadcast', 'Rowwise', 'I-ring' )
174 CALL pb_topset( ictxt, 'Broadcast', 'Columnwise', ' ' )
175
176 CALL descset( descv, m+iroff, desca( nb_ ), desca( mb_ ),
177 $ desca( nb_ ), iarow, iacol, ictxt,
max( 1, mp ) )
178
179
180
181 iv = 1 + m - k + iroff
182 jb = jn - ja - n + k + 1
183
184
185
186 CALL pzlarft(
'Backward',
'Columnwise', m-n+jn-ja+1, jb, a, ia,
187 $ ja+n-k, desca, tau, work( ipt ), work( ipw ) )
188
189
190
191 CALL pzlacpy(
'All', m-n+jn-ja+1, jb, a, ia, ja+n-k, desca,
192 $ work( ipv ), iroff+1, 1, descv )
193 CALL pzlaset(
'Lower', jb, jb, zero, one, work( ipv ), iv,
194 $ 1, descv )
195
196
197
198
199 CALL pzlaset(
'All', m-k, jb, zero, zero, a, ia, ja+n-k,
200 $ desca )
201 CALL pzlaset(
'Upper', jb, jb-1, zero, zero, a, ia+m-k,
202 $ ja+n-k+1, desca )
203
204
205
206 CALL pzlarfb(
'Left',
'No transpose',
'Backward',
'Columnwise',
207 $ m-n+jn-ja+1, jn-ja+1, jb, work( ipv ), iroff+1, 1,
208 $ descv, work( ipt ), a, ia, ja, desca, work( ipw ) )
209
210 descv( csrc_ ) = mod( descv( csrc_ ) + 1, npcol )
211
212
213
214 DO 10 j = jn+1, ja+n-1, desca( nb_ )
215 jb =
min( ja+n-j, desca( nb_ ) )
216 iv = 1 + m - n + j - ja + iroff
217
218
219
220 CALL pzlarft(
'Backward',
'Columnwise', m-n+j+jb-ja, jb, a, ia,
221 $ j, desca, tau, work( ipt ), work( ipw ) )
222
223
224
225 CALL pzlacpy(
'All', m-n+j+jb-ja, jb, a, ia, j, desca,
226 $ work( ipv ), iroff+1, 1, descv )
227 CALL pzlaset(
'Lower', jb, jb, zero, one, work( ipv ), iv,
228 $ 1, descv )
229
230
231
232
233 CALL pzlaset(
'All', m-n+j-ja, jb, zero, zero, a, ia, j,
234 $ desca )
235 CALL pzlaset(
'Upper', jb, jb-1, zero, zero, a, ia+m-n+j-ja,
236 $ j+1, desca )
237
238
239
240 CALL pzlarfb(
'Left',
'No transpose',
'Backward',
'Columnwise',
241 $ m-n+j+jb-ja, j+jb-ja, jb, work( ipv ), iroff+1,
242 $ 1, descv, work( ipt ), a, ia, ja, desca,
243 $ work( ipw ) )
244
245 descv( csrc_ ) = mod( descv( csrc_ ) + 1, npcol )
246
247 10 CONTINUE
248
249 CALL pb_topset( ictxt, 'Broadcast', 'Rowwise', rowbtop )
250 CALL pb_topset( ictxt, 'Broadcast', 'Columnwise', colbtop )
251
252 RETURN
253
254
255
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)