2
3
4
5
6
7
8
9 CHARACTER UPLO
10 INTEGER IA, JA, N
11
12
13 INTEGER DESCA( * )
14 COMPLEX*16 A( * )
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
128 parameter( one = ( 1.0d+0, 0.0d+0 ) )
129
130
131 INTEGER IACOL, IAROW, ICURR, IDIAG, IIA, IOFFA, JJA,
132 $ LDA, MYCOL, MYROW, NA, NPCOL, NPROW
133 DOUBLE PRECISION AII
134 COMPLEX*16 DOTC
135
136
137 EXTERNAL blacs_gridinfo,
infog2l, zdscal, zgemv,
139
140
141 LOGICAL LSAME
143
144
145 INTRINSIC dcmplx, dble
146
147
148
149
150
151 IF( n.EQ.0 )
152 $ RETURN
153
154
155
156 CALL blacs_gridinfo( desca( ctxt_ ), nprow, npcol, myrow, mycol )
157 CALL infog2l( ia, ja, desca, nprow, npcol, myrow, mycol, iia, jja,
158 $ iarow, iacol )
159
160 IF( myrow.EQ.iarow .AND. mycol.EQ.iacol ) THEN
161
162 lda = desca( lld_ )
163 idiag = iia + ( jja - 1 ) * lda
164 ioffa = idiag
165
166 IF(
lsame( uplo,
'U' ) )
THEN
167
168
169
170 DO 10 na = n-1, 1, -1
171 aii = a( idiag )
172 icurr = idiag + lda
173 CALL zzdotc( na, dotc, a( icurr ), lda, a( icurr ), lda )
174 a( idiag ) = aii*aii + dble( dotc )
175 CALL zlacgv( na, a( icurr ), lda )
176 CALL zgemv( 'No transpose', n-na-1, na, one,
177 $ a( ioffa+lda ), lda, a( icurr ), lda,
178 $ dcmplx( aii ), a( ioffa ), 1 )
179 CALL zlacgv( na, a( icurr ), lda )
180 idiag = idiag + lda + 1
181 ioffa = ioffa + lda
182 10 CONTINUE
183 aii = a( idiag )
184 CALL zdscal( n, aii, a( ioffa ), 1 )
185
186 ELSE
187
188
189
190 DO 20 na = 1, n-1
191 aii = a( idiag )
192 icurr = idiag + 1
193 a( idiag ) = aii*aii + dble( zdotc( n-na, a( icurr ), 1,
194 $ a( icurr ), 1 ) )
195 CALL zlacgv( na-1, a( ioffa ), lda )
196 CALL zgemv( 'Conjugate transpose', n-na, na-1, one,
197 $ a( ioffa+1 ), lda, a( icurr ), 1,
198 $ dcmplx( aii ), a( ioffa ), lda )
199 CALL zlacgv( na-1, a( ioffa ), lda )
200 idiag = idiag + lda + 1
201 ioffa = ioffa + 1
202 20 CONTINUE
203 aii = a( idiag )
204 CALL zdscal( n, aii, a( ioffa ), lda )
205
206 END IF
207
208 END IF
209
210 RETURN
211
212
213
subroutine infog2l(grindx, gcindx, desc, nprow, npcol, myrow, mycol, lrindx, lcindx, rsrc, csrc)
subroutine zzdotc(n, dotc, x, incx, y, incy)