4 SUBROUTINE d6dot(T,A,B,N)
10 DOUBLE PRECISION :: T(9)
11 DOUBLE PRECISION :: A(9,*)
12 DOUBLE PRECISION :: B(9,*)
25 t(1) = t(1) + a(1,jj)*b(1,jj) + a(4,jj)*b(4,jj) + a(7,jj)*b(7,jj)
26 t(2) = t(2) + a(2,jj)*b(1,jj) + a(5,jj)*b(4,jj) + a(8,jj)*b(7,jj)
27 t(3) = t(3) + a(3,jj)*b(1,jj) + a(6,jj)*b(4,jj) + a(9,jj)*b(7,jj)
28 t(4) = t(4) + a(1,jj)*b(2,jj) + a(4,jj)*b(5,jj) + a(7,jj)*b(8,jj)
29 t(5) = t(5) + a(2,jj)*b(2,jj) + a(5,jj)*b(5,jj) + a(8,jj)*b(8,jj)
30 t(6) = t(6) + a(3,jj)*b(2,jj) + a(6,jj)*b(5,jj) + a(9,jj)*b(8,jj)
31 t(7) = t(7) + a(1,jj)*b(3,jj) + a(4,jj)*b(6,jj) + a(7,jj)*b(9,jj)
32 t(8) = t(8) + a(2,jj)*b(3,jj) + a(5,jj)*b(6,jj) + a(8,jj)*b(9,jj)
33 t(9) = t(9) + a(3,jj)*b(3,jj) + a(6,jj)*b(6,jj) + a(9,jj)*b(9,jj)
45 DOUBLE PRECISION :: T(6)
46 DOUBLE PRECISION :: A(9,*)
47 DOUBLE PRECISION :: B(9,*)
61 t(1) = t(1) + a(1,jj)*b(1,jj) + a(4,jj)*b(4,jj) + a(7,jj)*b(7,jj)
62 t(2) = t(2) + a(2,jj)*b(1,jj) + a(5,jj)*b(4,jj) + a(8,jj)*b(7,jj)
63 t(3) = t(3) + a(2,jj)*b(2,jj) + a(5,jj)*b(5,jj) + a(8,jj)*b(8,jj)
64 t(4) = t(4) + a(3,jj)*b(1,jj) + a(6,jj)*b(4,jj) + a(9,jj)*b(7,jj)
65 t(5) = t(5) + a(3,jj)*b(2,jj) + a(6,jj)*b(5,jj) + a(9,jj)*b(8,jj)
66 t(6) = t(6) + a(3,jj)*b(3,jj) + a(6,jj)*b(6,jj) + a(9,jj)*b(9,jj)
72 SUBROUTINE d6sdot(Wi,A,B,N)
77 DOUBLE PRECISION :: Wi(3)
78 DOUBLE PRECISION :: A(3,*)
79 DOUBLE PRECISION :: B(9,*)
89 wi(1) = wi(1) - a(1,jj)*b(1,jj) - a(2,jj)*b(4,jj) - a(3,jj)*b(7,jj)
90 wi(2) = wi(2) - a(1,jj)*b(2,jj) - a(2,jj)*b(5,jj) - a(3,jj)*b(8,jj)
91 wi(3) = wi(3) - a(1,jj)*b(3,jj) - a(2,jj)*b(6,jj) - a(3,jj)*b(9,jj)
97 SUBROUTINE idntty(Neqns,Invp,Iperm)
105 COMMON /debug / idbg1
108 DO WHILE ( i<=neqns )
109 WRITE (6,*)
'invp(', i,
')'
111 IF ( invp(i)==0 )
THEN
117 ELSEIF ( invp(i)<0 )
THEN
118 READ (11,*) (invp(i),i=1,neqns)
131 99999
END SUBROUTINE idntty
135 SUBROUTINE nusol6(Xlnzr,Colno,Dsln,Zln,Diag,Iperm,B,Wk,Neqns,Nstop)
152 DOUBLE PRECISION :: Zln(36,*)
153 DOUBLE PRECISION :: Diag(21,*)
154 DOUBLE PRECISION :: B(6,*)
155 DOUBLE PRECISION :: Wk(6,*)
156 DOUBLE PRECISION :: Dsln(36,*)
159 wk(1,i) = b(1,iperm(i))
160 wk(2,i) = b(2,iperm(i))
161 wk(3,i) = b(3,iperm(i))
162 wk(4,i) = b(4,iperm(i))
163 wk(5,i) = b(5,iperm(i))
164 wk(6,i) = b(6,iperm(i))
170 IF ( ke>=ks )
CALL s6pdot(wk(1,i),wk,zln,colno,ks,ke)
173 CALL dxsdot(6,wk(1,i),wk(1,nstop),dsln(1,joc),i-nstop)
174 joc = joc + i - nstop
178 wk(2,i) = wk(2,i) - wk(1,i)*diag(2,i)
179 wk(3,i) = wk(3,i) - wk(1,i)*diag(4,i) - wk(2,i)*diag(5,i)
180 wk(4,i) = wk(4,i) - wk(1,i)*diag(7,i) - wk(2,i)*diag(8,i) - wk(3,i)*diag(9,i)
181 wk(5,i) = wk(5,i) - wk(1,i)*diag(11,i) - wk(2,i)*diag(12,i) - wk(3,i)*diag(13,i) - wk(4,i)*diag(14,i)
182 wk(6,i) = wk(6,i) - wk(1,i)*diag(16,i) - wk(2,i)*diag(17,i)&
183 - wk(3,i)*diag(18,i) - wk(4,i)*diag(19,i) - wk(6,i)&
185 wk(1,i) = wk(1,i)*diag(1,i)
186 wk(2,i) = wk(2,i)*diag(3,i)
187 wk(3,i) = wk(3,i)*diag(6,i)
188 wk(4,i) = wk(4,i)*diag(10,i)
189 wk(5,i) = wk(5,i)*diag(15,i)
190 wk(6,i) = wk(6,i)*diag(21,i)
191 wk(5,i) = wk(5,i) - wk(6,i)*diag(20,i)
192 wk(4,i) = wk(4,i) - wk(6,i)*diag(19,i) - wk(5,i)*diag(14,i)
193 wk(3,i) = wk(3,i) - wk(6,i)*diag(18,i) - wk(5,i)*diag(13,i) - wk(4,i)*diag(9,i)
194 wk(2,i) = wk(2,i) - wk(6,i)*diag(17,i) - wk(5,i)*diag(12,i) - wk(4,i)*diag(8,i) - wk(3,i)*diag(5,i)
195 wk(1,i) = wk(1,i) - wk(6,i)*diag(16,i) - wk(5,i)*diag(11,i)&
196 - wk(4,i)*diag(7,i) - wk(3,i)*diag(4,i) - wk(2,i)&
202 DO j = i - 1, nstop, -1
204 wk(1,j) = wk(1,j) - wk(1,i)*dsln(1,joc) - wk(2,i)&
205 *dsln(2,joc) - wk(3,i)*dsln(3,joc) - wk(4,i)&
206 *dsln(4,joc) - wk(5,i)*dsln(5,joc) - wk(6,i)&
208 wk(2,j) = wk(2,j) - wk(1,i)*dsln(7,joc) - wk(2,i)&
209 *dsln(8,joc) - wk(3,i)*dsln(9,joc) - wk(4,i)&
210 *dsln(10,joc) - wk(5,i)*dsln(11,joc) - wk(6,i)&
212 wk(3,j) = wk(3,j) - wk(1,i)*dsln(13,joc) - wk(2,i)&
213 *dsln(14,joc) - wk(3,i)*dsln(15,joc) - wk(4,i)&
214 *dsln(16,joc) - wk(5,i)*dsln(17,joc) - wk(6,i)&
216 wk(4,j) = wk(4,j) - wk(1,i)*dsln(19,joc) - wk(2,i)&
217 *dsln(20,joc) - wk(3,i)*dsln(21,joc) - wk(4,i)&
218 *dsln(22,joc) - wk(5,i)*dsln(23,joc) - wk(6,i)&
220 wk(5,j) = wk(5,j) - wk(1,i)*dsln(25,joc) - wk(2,i)&
221 *dsln(26,joc) - wk(3,i)*dsln(27,joc) - wk(4,i)&
222 *dsln(28,joc) - wk(5,i)*dsln(29,joc) - wk(6,i)&
224 wk(6,j) = wk(6,j) - wk(1,i)*dsln(31,joc) - wk(2,i)&
225 *dsln(32,joc) - wk(3,i)*dsln(33,joc) - wk(4,i)&
226 *dsln(34,joc) - wk(5,i)*dsln(35,joc) - wk(6,i)&
235 wk(1,j) = wk(1,j) - wk(1,i)*zln(1,joc) - wk(2,i)&
236 *zln(2,joc) - wk(3,i)*zln(3,joc) - wk(4,i)&
237 *zln(4,joc) - wk(5,i)*zln(5,joc) - wk(6,i)&
239 wk(2,j) = wk(2,j) - wk(1,i)*zln(7,joc) - wk(2,i)&
240 *zln(8,joc) - wk(3,i)*zln(9,joc) - wk(4,i)&
241 *zln(10,joc) - wk(5,i)*zln(11,joc) - wk(6,i)&
243 wk(3,j) = wk(3,j) - wk(1,i)*zln(13,joc) - wk(2,i)&
244 *zln(14,joc) - wk(3,i)*zln(15,joc) - wk(4,i)&
245 *zln(16,joc) - wk(5,i)*zln(17,joc) - wk(6,i)&
247 wk(4,j) = wk(4,j) - wk(1,i)*zln(19,joc) - wk(2,i)&
248 *zln(20,joc) - wk(3,i)*zln(21,joc) - wk(4,i)&
249 *zln(22,joc) - wk(5,i)*zln(23,joc) - wk(6,i)&
251 wk(5,j) = wk(5,j) - wk(1,i)*zln(25,joc) - wk(2,i)&
252 *zln(26,joc) - wk(3,i)*zln(27,joc) - wk(4,i)&
253 *zln(28,joc) - wk(5,i)*zln(29,joc) - wk(6,i)&
255 wk(6,j) = wk(6,j) - wk(1,i)*zln(31,joc) - wk(2,i)&
256 *zln(32,joc) - wk(3,i)*zln(33,joc) - wk(4,i)&
257 *zln(34,joc) - wk(5,i)*zln(35,joc) - wk(6,i)&
264 b(1,iperm(i)) = wk(1,i)
265 b(2,iperm(i)) = wk(2,i)
266 b(3,iperm(i)) = wk(3,i)
267 b(4,iperm(i)) = wk(4,i)
268 b(5,iperm(i)) = wk(5,i)
269 b(6,iperm(i)) = wk(6,i)
282 WRITE (6,99001) (ip(i),i=1,n)
283 99001
FORMAT (10(2x,i4))
289 DOUBLE PRECISION :: A
300 SUBROUTINE verif0(Neqns,Ndeg,Nttbr,Irow,Jcol,Val,Rhs,X)
303 DOUBLE PRECISION :: err
304 DOUBLE PRECISION :: rel
305 DOUBLE PRECISION :: Rhs
306 DOUBLE PRECISION :: Val
307 DOUBLE PRECISION :: X
318 dimension irow(*), jcol(*), val(ndeg,ndeg,*), rhs(ndeg,*), x(ndeg,*)
327 rel = rel + dabs(rhs(l,i))
335 rhs(l,i) = rhs(l,i) - val(l,m,k)*x(m,j)
336 IF ( i/=j ) rhs(l,j) = rhs(l,j) - val(m,l,k)*x(m,i)
343 err = err + dabs(rhs(l,i))
346 WRITE (6,99001) err, rel, err/rel
349 99001
FORMAT (
' ***verification***(symmetric)'/&
351 1pd20.10/
'norm(b) = ',&
352 1pd20.10/
'norm(Ax-b)/norm(b) = ',1pd20.10)
358 SUBROUTINE v6prod(Zln,Diag,Zz,N)
361 DOUBLE PRECISION :: Diag
362 DOUBLE PRECISION :: Zln
363 DOUBLE PRECISION :: Zz
366 dimension zln(9,n), diag(6,n), zz(9,n)
368 zz(4,i) = zln(4,i) - zln(1,i)*diag(2,i)
369 zz(7,i) = zln(7,i) - zln(1,i)*diag(4,i) - zz(4,i)*diag(5,i)
370 zz(1,i) = zln(1,i)*diag(1,i)
371 zz(4,i) = zz(4,i)*diag(3,i)
372 zz(7,i) = zz(7,i)*diag(6,i)
373 zz(4,i) = zz(4,i) - zz(7,i)*diag(5,i)
374 zz(1,i) = zz(1,i) - zz(4,i)*diag(2,i) - zz(7,i)*diag(4,i)
376 zz(5,i) = zln(5,i) - zln(2,i)*diag(2,i)
377 zz(8,i) = zln(8,i) - zln(2,i)*diag(4,i) - zz(5,i)*diag(5,i)
378 zz(2,i) = zln(2,i)*diag(1,i)
379 zz(5,i) = zz(5,i)*diag(3,i)
380 zz(8,i) = zz(8,i)*diag(6,i)
381 zz(5,i) = zz(5,i) - zz(8,i)*diag(5,i)
382 zz(2,i) = zz(2,i) - zz(5,i)*diag(2,i) - zz(8,i)*diag(4,i)
384 zz(6,i) = zln(6,i) - zln(3,i)*diag(2,i)
385 zz(9,i) = zln(9,i) - zln(3,i)*diag(4,i) - zz(6,i)*diag(5,i)
386 zz(3,i) = zln(3,i)*diag(1,i)
387 zz(6,i) = zz(6,i)*diag(3,i)
388 zz(9,i) = zz(9,i)*diag(6,i)
389 zz(6,i) = zz(6,i) - zz(9,i)*diag(5,i)
390 zz(3,i) = zz(3,i) - zz(6,i)*diag(2,i) - zz(9,i)*diag(4,i)