You can not select more than 25 topics Topics must start with a chinese character,a letter or number, can include dashes ('-') and can be up to 35 characters long.

ztrsvf.f 12 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361
  1. SUBROUTINE ZTRSVF ( UPLO, TRANS, DIAG, N, A, LDA, X, INCX )
  2. * .. Scalar Arguments ..
  3. INTEGER INCX, LDA, N
  4. CHARACTER*1 DIAG, TRANS, UPLO
  5. * .. Array Arguments ..
  6. COMPLEX*16 A( LDA, * ), X( * )
  7. * ..
  8. *
  9. * Purpose
  10. * =======
  11. *
  12. * ZTRSV solves one of the systems of equations
  13. *
  14. * A*x = b, or A'*x = b, or conjg( A' )*x = b,
  15. *
  16. * where b and x are n element vectors and A is an n by n unit, or
  17. * non-unit, upper or lower triangular matrix.
  18. *
  19. * No test for singularity or near-singularity is included in this
  20. * routine. Such tests must be performed before calling this routine.
  21. *
  22. * Parameters
  23. * ==========
  24. *
  25. * UPLO - CHARACTER*1.
  26. * On entry, UPLO specifies whether the matrix is an upper or
  27. * lower triangular matrix as follows:
  28. *
  29. * UPLO = 'U' or 'u' A is an upper triangular matrix.
  30. *
  31. * UPLO = 'L' or 'l' A is a lower triangular matrix.
  32. *
  33. * Unchanged on exit.
  34. *
  35. * TRANS - CHARACTER*1.
  36. * On entry, TRANS specifies the equations to be solved as
  37. * follows:
  38. *
  39. * TRANS = 'N' or 'n' A*x = b.
  40. *
  41. * TRANS = 'T' or 't' A'*x = b.
  42. *
  43. * TRANS = 'C' or 'c' conjg( A' )*x = b.
  44. *
  45. * Unchanged on exit.
  46. *
  47. * DIAG - CHARACTER*1.
  48. * On entry, DIAG specifies whether or not A is unit
  49. * triangular as follows:
  50. *
  51. * DIAG = 'U' or 'u' A is assumed to be unit triangular.
  52. *
  53. * DIAG = 'N' or 'n' A is not assumed to be unit
  54. * triangular.
  55. *
  56. * Unchanged on exit.
  57. *
  58. * N - INTEGER.
  59. * On entry, N specifies the order of the matrix A.
  60. * N must be at least zero.
  61. * Unchanged on exit.
  62. *
  63. * A - COMPLEX*16 array of DIMENSION ( LDA, n ).
  64. * Before entry with UPLO = 'U' or 'u', the leading n by n
  65. * upper triangular part of the array A must contain the upper
  66. * triangular matrix and the strictly lower triangular part of
  67. * A is not referenced.
  68. * Before entry with UPLO = 'L' or 'l', the leading n by n
  69. * lower triangular part of the array A must contain the lower
  70. * triangular matrix and the strictly upper triangular part of
  71. * A is not referenced.
  72. * Note that when DIAG = 'U' or 'u', the diagonal elements of
  73. * A are not referenced either, but are assumed to be unity.
  74. * Unchanged on exit.
  75. *
  76. * LDA - INTEGER.
  77. * On entry, LDA specifies the first dimension of A as declared
  78. * in the calling (sub) program. LDA must be at least
  79. * max( 1, n ).
  80. * Unchanged on exit.
  81. *
  82. * X - COMPLEX*16 array of dimension at least
  83. * ( 1 + ( n - 1 )*abs( INCX ) ).
  84. * Before entry, the incremented array X must contain the n
  85. * element right-hand side vector b. On exit, X is overwritten
  86. * with the solution vector x.
  87. *
  88. * INCX - INTEGER.
  89. * On entry, INCX specifies the increment for the elements of
  90. * X. INCX must not be zero.
  91. * Unchanged on exit.
  92. *
  93. *
  94. * Level 2 Blas routine.
  95. *
  96. * -- Written on 22-October-1986.
  97. * Jack Dongarra, Argonne National Lab.
  98. * Jeremy Du Croz, Nag Central Office.
  99. * Sven Hammarling, Nag Central Office.
  100. * Richard Hanson, Sandia National Labs.
  101. *
  102. *
  103. * .. Parameters ..
  104. COMPLEX*16 ZERO
  105. PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ) )
  106. * .. Local Scalars ..
  107. COMPLEX*16 TEMP
  108. INTEGER I, INFO, IX, J, JX, KX
  109. LOGICAL NOCONJ, NOUNIT
  110. * .. External Functions ..
  111. LOGICAL LSAME
  112. EXTERNAL LSAME
  113. * .. External Subroutines ..
  114. EXTERNAL XERBLA
  115. * .. Intrinsic Functions ..
  116. INTRINSIC DCONJG, MAX
  117. * ..
  118. * .. Executable Statements ..
  119. *
  120. * Test the input parameters.
  121. *
  122. INFO = 0
  123. IF ( .NOT.LSAME( UPLO , 'U' ).AND.
  124. $ .NOT.LSAME( UPLO , 'L' ) )THEN
  125. INFO = 1
  126. ELSE IF( .NOT.LSAME( TRANS, 'N' ).AND.
  127. $ .NOT.LSAME( TRANS, 'T' ).AND.
  128. $ .NOT.LSAME( TRANS, 'R' ).AND.
  129. $ .NOT.LSAME( TRANS, 'C' ) )THEN
  130. INFO = 2
  131. ELSE IF( .NOT.LSAME( DIAG , 'U' ).AND.
  132. $ .NOT.LSAME( DIAG , 'N' ) )THEN
  133. INFO = 3
  134. ELSE IF( N.LT.0 )THEN
  135. INFO = 4
  136. ELSE IF( LDA.LT.MAX( 1, N ) )THEN
  137. INFO = 6
  138. ELSE IF( INCX.EQ.0 )THEN
  139. INFO = 8
  140. END IF
  141. IF( INFO.NE.0 )THEN
  142. CALL XERBLA( 'ZTRSV ', INFO )
  143. RETURN
  144. END IF
  145. *
  146. * Quick return if possible.
  147. *
  148. IF( N.EQ.0 )
  149. $ RETURN
  150. *
  151. NOCONJ = LSAME( TRANS, 'N' ) .OR. LSAME( TRANS, 'T' )
  152. NOUNIT = LSAME( DIAG , 'N' )
  153. *
  154. * Set up the start point in X if the increment is not unity. This
  155. * will be ( N - 1 )*INCX too small for descending loops.
  156. *
  157. IF( INCX.LE.0 )THEN
  158. KX = 1 - ( N - 1 )*INCX
  159. ELSE IF( INCX.NE.1 )THEN
  160. KX = 1
  161. END IF
  162. *
  163. * Start the operations. In this version the elements of A are
  164. * accessed sequentially with one pass through A.
  165. *
  166. IF( LSAME( TRANS, 'N' ) .OR. LSAME( TRANS, 'R' ) ) THEN
  167. *
  168. * Form x := inv( A )*x.
  169. *
  170. IF( LSAME( UPLO, 'U' ) )THEN
  171. IF( INCX.EQ.1 )THEN
  172. DO 20, J = N, 1, -1
  173. IF( X( J ).NE.ZERO )THEN
  174. IF (NOCONJ) THEN
  175. IF( NOUNIT )
  176. $ X( J ) = X( J )/A( J, J )
  177. TEMP = X( J )
  178. DO 10, I = J - 1, 1, -1
  179. X( I ) = X( I ) - TEMP*A( I, J )
  180. 10 CONTINUE
  181. ELSE
  182. IF( NOUNIT )
  183. $ X( J ) = X( J )/DCONJG(A( J, J ))
  184. TEMP = X( J )
  185. DO 15, I = J - 1, 1, -1
  186. X( I ) = X( I ) - TEMP*DCONJG(A( I, J ))
  187. 15 CONTINUE
  188. ENDIF
  189. END IF
  190. 20 CONTINUE
  191. ELSE
  192. JX = KX + ( N - 1 )*INCX
  193. DO 40, J = N, 1, -1
  194. IF( X( JX ).NE.ZERO )THEN
  195. IF (NOCONJ) THEN
  196. IF( NOUNIT )
  197. $ X( JX ) = X( JX )/A( J, J )
  198. ELSE
  199. IF( NOUNIT )
  200. $ X( JX ) = X( JX )/DCONJG(A( J, J ))
  201. ENDIF
  202. TEMP = X( JX )
  203. IX = JX
  204. DO 30, I = J - 1, 1, -1
  205. IX = IX - INCX
  206. IF (NOCONJ) THEN
  207. X( IX ) = X( IX ) - TEMP*A( I, J )
  208. ELSE
  209. X( IX ) = X( IX ) - TEMP*DCONJG(A( I, J ))
  210. ENDIF
  211. 30 CONTINUE
  212. END IF
  213. JX = JX - INCX
  214. 40 CONTINUE
  215. END IF
  216. ELSE
  217. IF( INCX.EQ.1 )THEN
  218. DO 60, J = 1, N
  219. IF( X( J ).NE.ZERO )THEN
  220. IF (NOCONJ) THEN
  221. IF( NOUNIT )
  222. $ X( J ) = X( J )/A( J, J )
  223. TEMP = X( J )
  224. DO 50, I = J + 1, N
  225. X( I ) = X( I ) - TEMP*A( I, J )
  226. 50 CONTINUE
  227. ELSE
  228. IF( NOUNIT )
  229. $ X( J ) = X( J )/DCONJG(A( J, J ))
  230. TEMP = X( J )
  231. DO 55, I = J + 1, N
  232. X( I ) = X( I ) - TEMP*DCONJG(A( I, J ))
  233. 55 CONTINUE
  234. ENDIF
  235. END IF
  236. 60 CONTINUE
  237. ELSE
  238. JX = KX
  239. DO 80, J = 1, N
  240. IF( X( JX ).NE.ZERO )THEN
  241. IF (NOCONJ) THEN
  242. IF( NOUNIT )
  243. $ X( JX ) = X( JX )/A( J, J )
  244. ELSE
  245. IF( NOUNIT )
  246. $ X( JX ) = X( JX )/DCONJG(A( J, J ))
  247. ENDIF
  248. TEMP = X( JX )
  249. IX = JX
  250. DO 70, I = J + 1, N
  251. IX = IX + INCX
  252. IF (NOCONJ) THEN
  253. X( IX ) = X( IX ) - TEMP*A( I, J )
  254. ELSE
  255. X( IX ) = X( IX ) - TEMP*DCONJG(A( I, J ))
  256. ENDIF
  257. 70 CONTINUE
  258. END IF
  259. JX = JX + INCX
  260. 80 CONTINUE
  261. END IF
  262. END IF
  263. ELSE
  264. *
  265. * Form x := inv( A' )*x or x := inv( conjg( A' ) )*x.
  266. *
  267. IF( LSAME( UPLO, 'U' ) )THEN
  268. IF( INCX.EQ.1 )THEN
  269. DO 110, J = 1, N
  270. TEMP = X( J )
  271. IF( NOCONJ )THEN
  272. DO 90, I = 1, J - 1
  273. TEMP = TEMP - A( I, J )*X( I )
  274. 90 CONTINUE
  275. IF( NOUNIT )
  276. $ TEMP = TEMP/A( J, J )
  277. ELSE
  278. DO 100, I = 1, J - 1
  279. TEMP = TEMP - DCONJG( A( I, J ) )*X( I )
  280. 100 CONTINUE
  281. IF( NOUNIT )
  282. $ TEMP = TEMP/DCONJG( A( J, J ) )
  283. END IF
  284. X( J ) = TEMP
  285. 110 CONTINUE
  286. ELSE
  287. JX = KX
  288. DO 140, J = 1, N
  289. IX = KX
  290. TEMP = X( JX )
  291. IF( NOCONJ )THEN
  292. DO 120, I = 1, J - 1
  293. TEMP = TEMP - A( I, J )*X( IX )
  294. IX = IX + INCX
  295. 120 CONTINUE
  296. IF( NOUNIT )
  297. $ TEMP = TEMP/A( J, J )
  298. ELSE
  299. DO 130, I = 1, J - 1
  300. TEMP = TEMP - DCONJG( A( I, J ) )*X( IX )
  301. IX = IX + INCX
  302. 130 CONTINUE
  303. IF( NOUNIT )
  304. $ TEMP = TEMP/DCONJG( A( J, J ) )
  305. END IF
  306. X( JX ) = TEMP
  307. JX = JX + INCX
  308. 140 CONTINUE
  309. END IF
  310. ELSE
  311. IF( INCX.EQ.1 )THEN
  312. DO 170, J = N, 1, -1
  313. TEMP = X( J )
  314. IF( NOCONJ )THEN
  315. DO 150, I = N, J + 1, -1
  316. TEMP = TEMP - A( I, J )*X( I )
  317. 150 CONTINUE
  318. IF( NOUNIT )
  319. $ TEMP = TEMP/A( J, J )
  320. ELSE
  321. DO 160, I = N, J + 1, -1
  322. TEMP = TEMP - DCONJG( A( I, J ) )*X( I )
  323. 160 CONTINUE
  324. IF( NOUNIT )
  325. $ TEMP = TEMP/DCONJG( A( J, J ) )
  326. END IF
  327. X( J ) = TEMP
  328. 170 CONTINUE
  329. ELSE
  330. KX = KX + ( N - 1 )*INCX
  331. JX = KX
  332. DO 200, J = N, 1, -1
  333. IX = KX
  334. TEMP = X( JX )
  335. IF( NOCONJ )THEN
  336. DO 180, I = N, J + 1, -1
  337. TEMP = TEMP - A( I, J )*X( IX )
  338. IX = IX - INCX
  339. 180 CONTINUE
  340. IF( NOUNIT )
  341. $ TEMP = TEMP/A( J, J )
  342. ELSE
  343. DO 190, I = N, J + 1, -1
  344. TEMP = TEMP - DCONJG( A( I, J ) )*X( IX )
  345. IX = IX - INCX
  346. 190 CONTINUE
  347. IF( NOUNIT )
  348. $ TEMP = TEMP/DCONJG( A( J, J ) )
  349. END IF
  350. X( JX ) = TEMP
  351. JX = JX - INCX
  352. 200 CONTINUE
  353. END IF
  354. END IF
  355. END IF
  356. *
  357. RETURN
  358. *
  359. * End of ZTRSV .
  360. *
  361. END