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.

ctfsm.f 36 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026
  1. *> \brief \b CTFSM solves a matrix equation (one operand is a triangular matrix in RFP format).
  2. *
  3. * =========== DOCUMENTATION ===========
  4. *
  5. * Online html documentation available at
  6. * http://www.netlib.org/lapack/explore-html/
  7. *
  8. *> \htmlonly
  9. *> Download CTFSM + dependencies
  10. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ctfsm.f">
  11. *> [TGZ]</a>
  12. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ctfsm.f">
  13. *> [ZIP]</a>
  14. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ctfsm.f">
  15. *> [TXT]</a>
  16. *> \endhtmlonly
  17. *
  18. * Definition:
  19. * ===========
  20. *
  21. * SUBROUTINE CTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A,
  22. * B, LDB )
  23. *
  24. * .. Scalar Arguments ..
  25. * CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO
  26. * INTEGER LDB, M, N
  27. * COMPLEX ALPHA
  28. * ..
  29. * .. Array Arguments ..
  30. * COMPLEX A( 0: * ), B( 0: LDB-1, 0: * )
  31. * ..
  32. *
  33. *
  34. *> \par Purpose:
  35. * =============
  36. *>
  37. *> \verbatim
  38. *>
  39. *> Level 3 BLAS like routine for A in RFP Format.
  40. *>
  41. *> CTFSM solves the matrix equation
  42. *>
  43. *> op( A )*X = alpha*B or X*op( A ) = alpha*B
  44. *>
  45. *> where alpha is a scalar, X and B are m by n matrices, A is a unit, or
  46. *> non-unit, upper or lower triangular matrix and op( A ) is one of
  47. *>
  48. *> op( A ) = A or op( A ) = A**H.
  49. *>
  50. *> A is in Rectangular Full Packed (RFP) Format.
  51. *>
  52. *> The matrix X is overwritten on B.
  53. *> \endverbatim
  54. *
  55. * Arguments:
  56. * ==========
  57. *
  58. *> \param[in] TRANSR
  59. *> \verbatim
  60. *> TRANSR is CHARACTER*1
  61. *> = 'N': The Normal Form of RFP A is stored;
  62. *> = 'C': The Conjugate-transpose Form of RFP A is stored.
  63. *> \endverbatim
  64. *>
  65. *> \param[in] SIDE
  66. *> \verbatim
  67. *> SIDE is CHARACTER*1
  68. *> On entry, SIDE specifies whether op( A ) appears on the left
  69. *> or right of X as follows:
  70. *>
  71. *> SIDE = 'L' or 'l' op( A )*X = alpha*B.
  72. *>
  73. *> SIDE = 'R' or 'r' X*op( A ) = alpha*B.
  74. *>
  75. *> Unchanged on exit.
  76. *> \endverbatim
  77. *>
  78. *> \param[in] UPLO
  79. *> \verbatim
  80. *> UPLO is CHARACTER*1
  81. *> On entry, UPLO specifies whether the RFP matrix A came from
  82. *> an upper or lower triangular matrix as follows:
  83. *> UPLO = 'U' or 'u' RFP A came from an upper triangular matrix
  84. *> UPLO = 'L' or 'l' RFP A came from a lower triangular matrix
  85. *>
  86. *> Unchanged on exit.
  87. *> \endverbatim
  88. *>
  89. *> \param[in] TRANS
  90. *> \verbatim
  91. *> TRANS is CHARACTER*1
  92. *> On entry, TRANS specifies the form of op( A ) to be used
  93. *> in the matrix multiplication as follows:
  94. *>
  95. *> TRANS = 'N' or 'n' op( A ) = A.
  96. *>
  97. *> TRANS = 'C' or 'c' op( A ) = conjg( A' ).
  98. *>
  99. *> Unchanged on exit.
  100. *> \endverbatim
  101. *>
  102. *> \param[in] DIAG
  103. *> \verbatim
  104. *> DIAG is CHARACTER*1
  105. *> On entry, DIAG specifies whether or not RFP A is unit
  106. *> triangular as follows:
  107. *>
  108. *> DIAG = 'U' or 'u' A is assumed to be unit triangular.
  109. *>
  110. *> DIAG = 'N' or 'n' A is not assumed to be unit
  111. *> triangular.
  112. *>
  113. *> Unchanged on exit.
  114. *> \endverbatim
  115. *>
  116. *> \param[in] M
  117. *> \verbatim
  118. *> M is INTEGER
  119. *> On entry, M specifies the number of rows of B. M must be at
  120. *> least zero.
  121. *> Unchanged on exit.
  122. *> \endverbatim
  123. *>
  124. *> \param[in] N
  125. *> \verbatim
  126. *> N is INTEGER
  127. *> On entry, N specifies the number of columns of B. N must be
  128. *> at least zero.
  129. *> Unchanged on exit.
  130. *> \endverbatim
  131. *>
  132. *> \param[in] ALPHA
  133. *> \verbatim
  134. *> ALPHA is COMPLEX
  135. *> On entry, ALPHA specifies the scalar alpha. When alpha is
  136. *> zero then A is not referenced and B need not be set before
  137. *> entry.
  138. *> Unchanged on exit.
  139. *> \endverbatim
  140. *>
  141. *> \param[in] A
  142. *> \verbatim
  143. *> A is COMPLEX array, dimension (N*(N+1)/2)
  144. *> NT = N*(N+1)/2. On entry, the matrix A in RFP Format.
  145. *> RFP Format is described by TRANSR, UPLO and N as follows:
  146. *> If TRANSR='N' then RFP A is (0:N,0:K-1) when N is even;
  147. *> K=N/2. RFP A is (0:N-1,0:K) when N is odd; K=N/2. If
  148. *> TRANSR = 'C' then RFP is the Conjugate-transpose of RFP A as
  149. *> defined when TRANSR = 'N'. The contents of RFP A are defined
  150. *> by UPLO as follows: If UPLO = 'U' the RFP A contains the NT
  151. *> elements of upper packed A either in normal or
  152. *> conjugate-transpose Format. If UPLO = 'L' the RFP A contains
  153. *> the NT elements of lower packed A either in normal or
  154. *> conjugate-transpose Format. The LDA of RFP A is (N+1)/2 when
  155. *> TRANSR = 'C'. When TRANSR is 'N' the LDA is N+1 when N is
  156. *> even and is N when is odd.
  157. *> See the Note below for more details. Unchanged on exit.
  158. *> \endverbatim
  159. *>
  160. *> \param[in,out] B
  161. *> \verbatim
  162. *> B is COMPLEX array, dimension (LDB,N)
  163. *> Before entry, the leading m by n part of the array B must
  164. *> contain the right-hand side matrix B, and on exit is
  165. *> overwritten by the solution matrix X.
  166. *> \endverbatim
  167. *>
  168. *> \param[in] LDB
  169. *> \verbatim
  170. *> LDB is INTEGER
  171. *> On entry, LDB specifies the first dimension of B as declared
  172. *> in the calling (sub) program. LDB must be at least
  173. *> max( 1, m ).
  174. *> Unchanged on exit.
  175. *> \endverbatim
  176. *
  177. * Authors:
  178. * ========
  179. *
  180. *> \author Univ. of Tennessee
  181. *> \author Univ. of California Berkeley
  182. *> \author Univ. of Colorado Denver
  183. *> \author NAG Ltd.
  184. *
  185. *> \date December 2016
  186. *
  187. *> \ingroup complexOTHERcomputational
  188. *
  189. *> \par Further Details:
  190. * =====================
  191. *>
  192. *> \verbatim
  193. *>
  194. *> We first consider Standard Packed Format when N is even.
  195. *> We give an example where N = 6.
  196. *>
  197. *> AP is Upper AP is Lower
  198. *>
  199. *> 00 01 02 03 04 05 00
  200. *> 11 12 13 14 15 10 11
  201. *> 22 23 24 25 20 21 22
  202. *> 33 34 35 30 31 32 33
  203. *> 44 45 40 41 42 43 44
  204. *> 55 50 51 52 53 54 55
  205. *>
  206. *>
  207. *> Let TRANSR = 'N'. RFP holds AP as follows:
  208. *> For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
  209. *> three columns of AP upper. The lower triangle A(4:6,0:2) consists of
  210. *> conjugate-transpose of the first three columns of AP upper.
  211. *> For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
  212. *> three columns of AP lower. The upper triangle A(0:2,0:2) consists of
  213. *> conjugate-transpose of the last three columns of AP lower.
  214. *> To denote conjugate we place -- above the element. This covers the
  215. *> case N even and TRANSR = 'N'.
  216. *>
  217. *> RFP A RFP A
  218. *>
  219. *> -- -- --
  220. *> 03 04 05 33 43 53
  221. *> -- --
  222. *> 13 14 15 00 44 54
  223. *> --
  224. *> 23 24 25 10 11 55
  225. *>
  226. *> 33 34 35 20 21 22
  227. *> --
  228. *> 00 44 45 30 31 32
  229. *> -- --
  230. *> 01 11 55 40 41 42
  231. *> -- -- --
  232. *> 02 12 22 50 51 52
  233. *>
  234. *> Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
  235. *> transpose of RFP A above. One therefore gets:
  236. *>
  237. *>
  238. *> RFP A RFP A
  239. *>
  240. *> -- -- -- -- -- -- -- -- -- --
  241. *> 03 13 23 33 00 01 02 33 00 10 20 30 40 50
  242. *> -- -- -- -- -- -- -- -- -- --
  243. *> 04 14 24 34 44 11 12 43 44 11 21 31 41 51
  244. *> -- -- -- -- -- -- -- -- -- --
  245. *> 05 15 25 35 45 55 22 53 54 55 22 32 42 52
  246. *>
  247. *>
  248. *> We next consider Standard Packed Format when N is odd.
  249. *> We give an example where N = 5.
  250. *>
  251. *> AP is Upper AP is Lower
  252. *>
  253. *> 00 01 02 03 04 00
  254. *> 11 12 13 14 10 11
  255. *> 22 23 24 20 21 22
  256. *> 33 34 30 31 32 33
  257. *> 44 40 41 42 43 44
  258. *>
  259. *>
  260. *> Let TRANSR = 'N'. RFP holds AP as follows:
  261. *> For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
  262. *> three columns of AP upper. The lower triangle A(3:4,0:1) consists of
  263. *> conjugate-transpose of the first two columns of AP upper.
  264. *> For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
  265. *> three columns of AP lower. The upper triangle A(0:1,1:2) consists of
  266. *> conjugate-transpose of the last two columns of AP lower.
  267. *> To denote conjugate we place -- above the element. This covers the
  268. *> case N odd and TRANSR = 'N'.
  269. *>
  270. *> RFP A RFP A
  271. *>
  272. *> -- --
  273. *> 02 03 04 00 33 43
  274. *> --
  275. *> 12 13 14 10 11 44
  276. *>
  277. *> 22 23 24 20 21 22
  278. *> --
  279. *> 00 33 34 30 31 32
  280. *> -- --
  281. *> 01 11 44 40 41 42
  282. *>
  283. *> Now let TRANSR = 'C'. RFP A in both UPLO cases is just the conjugate-
  284. *> transpose of RFP A above. One therefore gets:
  285. *>
  286. *>
  287. *> RFP A RFP A
  288. *>
  289. *> -- -- -- -- -- -- -- -- --
  290. *> 02 12 22 00 01 00 10 20 30 40 50
  291. *> -- -- -- -- -- -- -- -- --
  292. *> 03 13 23 33 11 33 11 21 31 41 51
  293. *> -- -- -- -- -- -- -- -- --
  294. *> 04 14 24 34 44 43 44 22 32 42 52
  295. *> \endverbatim
  296. *>
  297. * =====================================================================
  298. SUBROUTINE CTFSM( TRANSR, SIDE, UPLO, TRANS, DIAG, M, N, ALPHA, A,
  299. $ B, LDB )
  300. *
  301. * -- LAPACK computational routine (version 3.7.0) --
  302. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  303. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  304. * December 2016
  305. *
  306. * .. Scalar Arguments ..
  307. CHARACTER TRANSR, DIAG, SIDE, TRANS, UPLO
  308. INTEGER LDB, M, N
  309. COMPLEX ALPHA
  310. * ..
  311. * .. Array Arguments ..
  312. COMPLEX A( 0: * ), B( 0: LDB-1, 0: * )
  313. * ..
  314. *
  315. * =====================================================================
  316. * ..
  317. * .. Parameters ..
  318. COMPLEX CONE, CZERO
  319. PARAMETER ( CONE = ( 1.0E+0, 0.0E+0 ),
  320. $ CZERO = ( 0.0E+0, 0.0E+0 ) )
  321. * ..
  322. * .. Local Scalars ..
  323. LOGICAL LOWER, LSIDE, MISODD, NISODD, NORMALTRANSR,
  324. $ NOTRANS
  325. INTEGER M1, M2, N1, N2, K, INFO, I, J
  326. * ..
  327. * .. External Functions ..
  328. LOGICAL LSAME
  329. EXTERNAL LSAME
  330. * ..
  331. * .. External Subroutines ..
  332. EXTERNAL XERBLA, CGEMM, CTRSM
  333. * ..
  334. * .. Intrinsic Functions ..
  335. INTRINSIC MAX, MOD
  336. * ..
  337. * .. Executable Statements ..
  338. *
  339. * Test the input parameters.
  340. *
  341. INFO = 0
  342. NORMALTRANSR = LSAME( TRANSR, 'N' )
  343. LSIDE = LSAME( SIDE, 'L' )
  344. LOWER = LSAME( UPLO, 'L' )
  345. NOTRANS = LSAME( TRANS, 'N' )
  346. IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'C' ) ) THEN
  347. INFO = -1
  348. ELSE IF( .NOT.LSIDE .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
  349. INFO = -2
  350. ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
  351. INFO = -3
  352. ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'C' ) ) THEN
  353. INFO = -4
  354. ELSE IF( .NOT.LSAME( DIAG, 'N' ) .AND. .NOT.LSAME( DIAG, 'U' ) )
  355. $ THEN
  356. INFO = -5
  357. ELSE IF( M.LT.0 ) THEN
  358. INFO = -6
  359. ELSE IF( N.LT.0 ) THEN
  360. INFO = -7
  361. ELSE IF( LDB.LT.MAX( 1, M ) ) THEN
  362. INFO = -11
  363. END IF
  364. IF( INFO.NE.0 ) THEN
  365. CALL XERBLA( 'CTFSM ', -INFO )
  366. RETURN
  367. END IF
  368. *
  369. * Quick return when ( (N.EQ.0).OR.(M.EQ.0) )
  370. *
  371. IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) )
  372. $ RETURN
  373. *
  374. * Quick return when ALPHA.EQ.(0E+0,0E+0)
  375. *
  376. IF( ALPHA.EQ.CZERO ) THEN
  377. DO 20 J = 0, N - 1
  378. DO 10 I = 0, M - 1
  379. B( I, J ) = CZERO
  380. 10 CONTINUE
  381. 20 CONTINUE
  382. RETURN
  383. END IF
  384. *
  385. IF( LSIDE ) THEN
  386. *
  387. * SIDE = 'L'
  388. *
  389. * A is M-by-M.
  390. * If M is odd, set NISODD = .TRUE., and M1 and M2.
  391. * If M is even, NISODD = .FALSE., and M.
  392. *
  393. IF( MOD( M, 2 ).EQ.0 ) THEN
  394. MISODD = .FALSE.
  395. K = M / 2
  396. ELSE
  397. MISODD = .TRUE.
  398. IF( LOWER ) THEN
  399. M2 = M / 2
  400. M1 = M - M2
  401. ELSE
  402. M1 = M / 2
  403. M2 = M - M1
  404. END IF
  405. END IF
  406. *
  407. IF( MISODD ) THEN
  408. *
  409. * SIDE = 'L' and N is odd
  410. *
  411. IF( NORMALTRANSR ) THEN
  412. *
  413. * SIDE = 'L', N is odd, and TRANSR = 'N'
  414. *
  415. IF( LOWER ) THEN
  416. *
  417. * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'L'
  418. *
  419. IF( NOTRANS ) THEN
  420. *
  421. * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
  422. * TRANS = 'N'
  423. *
  424. IF( M.EQ.1 ) THEN
  425. CALL CTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
  426. $ A, M, B, LDB )
  427. ELSE
  428. CALL CTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
  429. $ A( 0 ), M, B, LDB )
  430. CALL CGEMM( 'N', 'N', M2, N, M1, -CONE, A( M1 ),
  431. $ M, B, LDB, ALPHA, B( M1, 0 ), LDB )
  432. CALL CTRSM( 'L', 'U', 'C', DIAG, M2, N, CONE,
  433. $ A( M ), M, B( M1, 0 ), LDB )
  434. END IF
  435. *
  436. ELSE
  437. *
  438. * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and
  439. * TRANS = 'C'
  440. *
  441. IF( M.EQ.1 ) THEN
  442. CALL CTRSM( 'L', 'L', 'C', DIAG, M1, N, ALPHA,
  443. $ A( 0 ), M, B, LDB )
  444. ELSE
  445. CALL CTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
  446. $ A( M ), M, B( M1, 0 ), LDB )
  447. CALL CGEMM( 'C', 'N', M1, N, M2, -CONE, A( M1 ),
  448. $ M, B( M1, 0 ), LDB, ALPHA, B, LDB )
  449. CALL CTRSM( 'L', 'L', 'C', DIAG, M1, N, CONE,
  450. $ A( 0 ), M, B, LDB )
  451. END IF
  452. *
  453. END IF
  454. *
  455. ELSE
  456. *
  457. * SIDE ='L', N is odd, TRANSR = 'N', and UPLO = 'U'
  458. *
  459. IF( .NOT.NOTRANS ) THEN
  460. *
  461. * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
  462. * TRANS = 'N'
  463. *
  464. CALL CTRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA,
  465. $ A( M2 ), M, B, LDB )
  466. CALL CGEMM( 'C', 'N', M2, N, M1, -CONE, A( 0 ), M,
  467. $ B, LDB, ALPHA, B( M1, 0 ), LDB )
  468. CALL CTRSM( 'L', 'U', 'C', DIAG, M2, N, CONE,
  469. $ A( M1 ), M, B( M1, 0 ), LDB )
  470. *
  471. ELSE
  472. *
  473. * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'U', and
  474. * TRANS = 'C'
  475. *
  476. CALL CTRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA,
  477. $ A( M1 ), M, B( M1, 0 ), LDB )
  478. CALL CGEMM( 'N', 'N', M1, N, M2, -CONE, A( 0 ), M,
  479. $ B( M1, 0 ), LDB, ALPHA, B, LDB )
  480. CALL CTRSM( 'L', 'L', 'C', DIAG, M1, N, CONE,
  481. $ A( M2 ), M, B, LDB )
  482. *
  483. END IF
  484. *
  485. END IF
  486. *
  487. ELSE
  488. *
  489. * SIDE = 'L', N is odd, and TRANSR = 'C'
  490. *
  491. IF( LOWER ) THEN
  492. *
  493. * SIDE ='L', N is odd, TRANSR = 'C', and UPLO = 'L'
  494. *
  495. IF( NOTRANS ) THEN
  496. *
  497. * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'L', and
  498. * TRANS = 'N'
  499. *
  500. IF( M.EQ.1 ) THEN
  501. CALL CTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
  502. $ A( 0 ), M1, B, LDB )
  503. ELSE
  504. CALL CTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
  505. $ A( 0 ), M1, B, LDB )
  506. CALL CGEMM( 'C', 'N', M2, N, M1, -CONE,
  507. $ A( M1*M1 ), M1, B, LDB, ALPHA,
  508. $ B( M1, 0 ), LDB )
  509. CALL CTRSM( 'L', 'L', 'N', DIAG, M2, N, CONE,
  510. $ A( 1 ), M1, B( M1, 0 ), LDB )
  511. END IF
  512. *
  513. ELSE
  514. *
  515. * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'L', and
  516. * TRANS = 'C'
  517. *
  518. IF( M.EQ.1 ) THEN
  519. CALL CTRSM( 'L', 'U', 'N', DIAG, M1, N, ALPHA,
  520. $ A( 0 ), M1, B, LDB )
  521. ELSE
  522. CALL CTRSM( 'L', 'L', 'C', DIAG, M2, N, ALPHA,
  523. $ A( 1 ), M1, B( M1, 0 ), LDB )
  524. CALL CGEMM( 'N', 'N', M1, N, M2, -CONE,
  525. $ A( M1*M1 ), M1, B( M1, 0 ), LDB,
  526. $ ALPHA, B, LDB )
  527. CALL CTRSM( 'L', 'U', 'N', DIAG, M1, N, CONE,
  528. $ A( 0 ), M1, B, LDB )
  529. END IF
  530. *
  531. END IF
  532. *
  533. ELSE
  534. *
  535. * SIDE ='L', N is odd, TRANSR = 'C', and UPLO = 'U'
  536. *
  537. IF( .NOT.NOTRANS ) THEN
  538. *
  539. * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'U', and
  540. * TRANS = 'N'
  541. *
  542. CALL CTRSM( 'L', 'U', 'C', DIAG, M1, N, ALPHA,
  543. $ A( M2*M2 ), M2, B, LDB )
  544. CALL CGEMM( 'N', 'N', M2, N, M1, -CONE, A( 0 ), M2,
  545. $ B, LDB, ALPHA, B( M1, 0 ), LDB )
  546. CALL CTRSM( 'L', 'L', 'N', DIAG, M2, N, CONE,
  547. $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
  548. *
  549. ELSE
  550. *
  551. * SIDE ='L', N is odd, TRANSR = 'C', UPLO = 'U', and
  552. * TRANS = 'C'
  553. *
  554. CALL CTRSM( 'L', 'L', 'C', DIAG, M2, N, ALPHA,
  555. $ A( M1*M2 ), M2, B( M1, 0 ), LDB )
  556. CALL CGEMM( 'C', 'N', M1, N, M2, -CONE, A( 0 ), M2,
  557. $ B( M1, 0 ), LDB, ALPHA, B, LDB )
  558. CALL CTRSM( 'L', 'U', 'N', DIAG, M1, N, CONE,
  559. $ A( M2*M2 ), M2, B, LDB )
  560. *
  561. END IF
  562. *
  563. END IF
  564. *
  565. END IF
  566. *
  567. ELSE
  568. *
  569. * SIDE = 'L' and N is even
  570. *
  571. IF( NORMALTRANSR ) THEN
  572. *
  573. * SIDE = 'L', N is even, and TRANSR = 'N'
  574. *
  575. IF( LOWER ) THEN
  576. *
  577. * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'L'
  578. *
  579. IF( NOTRANS ) THEN
  580. *
  581. * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
  582. * and TRANS = 'N'
  583. *
  584. CALL CTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
  585. $ A( 1 ), M+1, B, LDB )
  586. CALL CGEMM( 'N', 'N', K, N, K, -CONE, A( K+1 ),
  587. $ M+1, B, LDB, ALPHA, B( K, 0 ), LDB )
  588. CALL CTRSM( 'L', 'U', 'C', DIAG, K, N, CONE,
  589. $ A( 0 ), M+1, B( K, 0 ), LDB )
  590. *
  591. ELSE
  592. *
  593. * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'L',
  594. * and TRANS = 'C'
  595. *
  596. CALL CTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
  597. $ A( 0 ), M+1, B( K, 0 ), LDB )
  598. CALL CGEMM( 'C', 'N', K, N, K, -CONE, A( K+1 ),
  599. $ M+1, B( K, 0 ), LDB, ALPHA, B, LDB )
  600. CALL CTRSM( 'L', 'L', 'C', DIAG, K, N, CONE,
  601. $ A( 1 ), M+1, B, LDB )
  602. *
  603. END IF
  604. *
  605. ELSE
  606. *
  607. * SIDE ='L', N is even, TRANSR = 'N', and UPLO = 'U'
  608. *
  609. IF( .NOT.NOTRANS ) THEN
  610. *
  611. * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
  612. * and TRANS = 'N'
  613. *
  614. CALL CTRSM( 'L', 'L', 'N', DIAG, K, N, ALPHA,
  615. $ A( K+1 ), M+1, B, LDB )
  616. CALL CGEMM( 'C', 'N', K, N, K, -CONE, A( 0 ), M+1,
  617. $ B, LDB, ALPHA, B( K, 0 ), LDB )
  618. CALL CTRSM( 'L', 'U', 'C', DIAG, K, N, CONE,
  619. $ A( K ), M+1, B( K, 0 ), LDB )
  620. *
  621. ELSE
  622. *
  623. * SIDE ='L', N is even, TRANSR = 'N', UPLO = 'U',
  624. * and TRANS = 'C'
  625. CALL CTRSM( 'L', 'U', 'N', DIAG, K, N, ALPHA,
  626. $ A( K ), M+1, B( K, 0 ), LDB )
  627. CALL CGEMM( 'N', 'N', K, N, K, -CONE, A( 0 ), M+1,
  628. $ B( K, 0 ), LDB, ALPHA, B, LDB )
  629. CALL CTRSM( 'L', 'L', 'C', DIAG, K, N, CONE,
  630. $ A( K+1 ), M+1, B, LDB )
  631. *
  632. END IF
  633. *
  634. END IF
  635. *
  636. ELSE
  637. *
  638. * SIDE = 'L', N is even, and TRANSR = 'C'
  639. *
  640. IF( LOWER ) THEN
  641. *
  642. * SIDE ='L', N is even, TRANSR = 'C', and UPLO = 'L'
  643. *
  644. IF( NOTRANS ) THEN
  645. *
  646. * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'L',
  647. * and TRANS = 'N'
  648. *
  649. CALL CTRSM( 'L', 'U', 'C', DIAG, K, N, ALPHA,
  650. $ A( K ), K, B, LDB )
  651. CALL CGEMM( 'C', 'N', K, N, K, -CONE,
  652. $ A( K*( K+1 ) ), K, B, LDB, ALPHA,
  653. $ B( K, 0 ), LDB )
  654. CALL CTRSM( 'L', 'L', 'N', DIAG, K, N, CONE,
  655. $ A( 0 ), K, B( K, 0 ), LDB )
  656. *
  657. ELSE
  658. *
  659. * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'L',
  660. * and TRANS = 'C'
  661. *
  662. CALL CTRSM( 'L', 'L', 'C', DIAG, K, N, ALPHA,
  663. $ A( 0 ), K, B( K, 0 ), LDB )
  664. CALL CGEMM( 'N', 'N', K, N, K, -CONE,
  665. $ A( K*( K+1 ) ), K, B( K, 0 ), LDB,
  666. $ ALPHA, B, LDB )
  667. CALL CTRSM( 'L', 'U', 'N', DIAG, K, N, CONE,
  668. $ A( K ), K, B, LDB )
  669. *
  670. END IF
  671. *
  672. ELSE
  673. *
  674. * SIDE ='L', N is even, TRANSR = 'C', and UPLO = 'U'
  675. *
  676. IF( .NOT.NOTRANS ) THEN
  677. *
  678. * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'U',
  679. * and TRANS = 'N'
  680. *
  681. CALL CTRSM( 'L', 'U', 'C', DIAG, K, N, ALPHA,
  682. $ A( K*( K+1 ) ), K, B, LDB )
  683. CALL CGEMM( 'N', 'N', K, N, K, -CONE, A( 0 ), K, B,
  684. $ LDB, ALPHA, B( K, 0 ), LDB )
  685. CALL CTRSM( 'L', 'L', 'N', DIAG, K, N, CONE,
  686. $ A( K*K ), K, B( K, 0 ), LDB )
  687. *
  688. ELSE
  689. *
  690. * SIDE ='L', N is even, TRANSR = 'C', UPLO = 'U',
  691. * and TRANS = 'C'
  692. *
  693. CALL CTRSM( 'L', 'L', 'C', DIAG, K, N, ALPHA,
  694. $ A( K*K ), K, B( K, 0 ), LDB )
  695. CALL CGEMM( 'C', 'N', K, N, K, -CONE, A( 0 ), K,
  696. $ B( K, 0 ), LDB, ALPHA, B, LDB )
  697. CALL CTRSM( 'L', 'U', 'N', DIAG, K, N, CONE,
  698. $ A( K*( K+1 ) ), K, B, LDB )
  699. *
  700. END IF
  701. *
  702. END IF
  703. *
  704. END IF
  705. *
  706. END IF
  707. *
  708. ELSE
  709. *
  710. * SIDE = 'R'
  711. *
  712. * A is N-by-N.
  713. * If N is odd, set NISODD = .TRUE., and N1 and N2.
  714. * If N is even, NISODD = .FALSE., and K.
  715. *
  716. IF( MOD( N, 2 ).EQ.0 ) THEN
  717. NISODD = .FALSE.
  718. K = N / 2
  719. ELSE
  720. NISODD = .TRUE.
  721. IF( LOWER ) THEN
  722. N2 = N / 2
  723. N1 = N - N2
  724. ELSE
  725. N1 = N / 2
  726. N2 = N - N1
  727. END IF
  728. END IF
  729. *
  730. IF( NISODD ) THEN
  731. *
  732. * SIDE = 'R' and N is odd
  733. *
  734. IF( NORMALTRANSR ) THEN
  735. *
  736. * SIDE = 'R', N is odd, and TRANSR = 'N'
  737. *
  738. IF( LOWER ) THEN
  739. *
  740. * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'L'
  741. *
  742. IF( NOTRANS ) THEN
  743. *
  744. * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
  745. * TRANS = 'N'
  746. *
  747. CALL CTRSM( 'R', 'U', 'C', DIAG, M, N2, ALPHA,
  748. $ A( N ), N, B( 0, N1 ), LDB )
  749. CALL CGEMM( 'N', 'N', M, N1, N2, -CONE, B( 0, N1 ),
  750. $ LDB, A( N1 ), N, ALPHA, B( 0, 0 ),
  751. $ LDB )
  752. CALL CTRSM( 'R', 'L', 'N', DIAG, M, N1, CONE,
  753. $ A( 0 ), N, B( 0, 0 ), LDB )
  754. *
  755. ELSE
  756. *
  757. * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'L', and
  758. * TRANS = 'C'
  759. *
  760. CALL CTRSM( 'R', 'L', 'C', DIAG, M, N1, ALPHA,
  761. $ A( 0 ), N, B( 0, 0 ), LDB )
  762. CALL CGEMM( 'N', 'C', M, N2, N1, -CONE, B( 0, 0 ),
  763. $ LDB, A( N1 ), N, ALPHA, B( 0, N1 ),
  764. $ LDB )
  765. CALL CTRSM( 'R', 'U', 'N', DIAG, M, N2, CONE,
  766. $ A( N ), N, B( 0, N1 ), LDB )
  767. *
  768. END IF
  769. *
  770. ELSE
  771. *
  772. * SIDE ='R', N is odd, TRANSR = 'N', and UPLO = 'U'
  773. *
  774. IF( NOTRANS ) THEN
  775. *
  776. * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
  777. * TRANS = 'N'
  778. *
  779. CALL CTRSM( 'R', 'L', 'C', DIAG, M, N1, ALPHA,
  780. $ A( N2 ), N, B( 0, 0 ), LDB )
  781. CALL CGEMM( 'N', 'N', M, N2, N1, -CONE, B( 0, 0 ),
  782. $ LDB, A( 0 ), N, ALPHA, B( 0, N1 ),
  783. $ LDB )
  784. CALL CTRSM( 'R', 'U', 'N', DIAG, M, N2, CONE,
  785. $ A( N1 ), N, B( 0, N1 ), LDB )
  786. *
  787. ELSE
  788. *
  789. * SIDE ='R', N is odd, TRANSR = 'N', UPLO = 'U', and
  790. * TRANS = 'C'
  791. *
  792. CALL CTRSM( 'R', 'U', 'C', DIAG, M, N2, ALPHA,
  793. $ A( N1 ), N, B( 0, N1 ), LDB )
  794. CALL CGEMM( 'N', 'C', M, N1, N2, -CONE, B( 0, N1 ),
  795. $ LDB, A( 0 ), N, ALPHA, B( 0, 0 ), LDB )
  796. CALL CTRSM( 'R', 'L', 'N', DIAG, M, N1, CONE,
  797. $ A( N2 ), N, B( 0, 0 ), LDB )
  798. *
  799. END IF
  800. *
  801. END IF
  802. *
  803. ELSE
  804. *
  805. * SIDE = 'R', N is odd, and TRANSR = 'C'
  806. *
  807. IF( LOWER ) THEN
  808. *
  809. * SIDE ='R', N is odd, TRANSR = 'C', and UPLO = 'L'
  810. *
  811. IF( NOTRANS ) THEN
  812. *
  813. * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'L', and
  814. * TRANS = 'N'
  815. *
  816. CALL CTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
  817. $ A( 1 ), N1, B( 0, N1 ), LDB )
  818. CALL CGEMM( 'N', 'C', M, N1, N2, -CONE, B( 0, N1 ),
  819. $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, 0 ),
  820. $ LDB )
  821. CALL CTRSM( 'R', 'U', 'C', DIAG, M, N1, CONE,
  822. $ A( 0 ), N1, B( 0, 0 ), LDB )
  823. *
  824. ELSE
  825. *
  826. * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'L', and
  827. * TRANS = 'C'
  828. *
  829. CALL CTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
  830. $ A( 0 ), N1, B( 0, 0 ), LDB )
  831. CALL CGEMM( 'N', 'N', M, N2, N1, -CONE, B( 0, 0 ),
  832. $ LDB, A( N1*N1 ), N1, ALPHA, B( 0, N1 ),
  833. $ LDB )
  834. CALL CTRSM( 'R', 'L', 'C', DIAG, M, N2, CONE,
  835. $ A( 1 ), N1, B( 0, N1 ), LDB )
  836. *
  837. END IF
  838. *
  839. ELSE
  840. *
  841. * SIDE ='R', N is odd, TRANSR = 'C', and UPLO = 'U'
  842. *
  843. IF( NOTRANS ) THEN
  844. *
  845. * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'U', and
  846. * TRANS = 'N'
  847. *
  848. CALL CTRSM( 'R', 'U', 'N', DIAG, M, N1, ALPHA,
  849. $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
  850. CALL CGEMM( 'N', 'C', M, N2, N1, -CONE, B( 0, 0 ),
  851. $ LDB, A( 0 ), N2, ALPHA, B( 0, N1 ),
  852. $ LDB )
  853. CALL CTRSM( 'R', 'L', 'C', DIAG, M, N2, CONE,
  854. $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
  855. *
  856. ELSE
  857. *
  858. * SIDE ='R', N is odd, TRANSR = 'C', UPLO = 'U', and
  859. * TRANS = 'C'
  860. *
  861. CALL CTRSM( 'R', 'L', 'N', DIAG, M, N2, ALPHA,
  862. $ A( N1*N2 ), N2, B( 0, N1 ), LDB )
  863. CALL CGEMM( 'N', 'N', M, N1, N2, -CONE, B( 0, N1 ),
  864. $ LDB, A( 0 ), N2, ALPHA, B( 0, 0 ),
  865. $ LDB )
  866. CALL CTRSM( 'R', 'U', 'C', DIAG, M, N1, CONE,
  867. $ A( N2*N2 ), N2, B( 0, 0 ), LDB )
  868. *
  869. END IF
  870. *
  871. END IF
  872. *
  873. END IF
  874. *
  875. ELSE
  876. *
  877. * SIDE = 'R' and N is even
  878. *
  879. IF( NORMALTRANSR ) THEN
  880. *
  881. * SIDE = 'R', N is even, and TRANSR = 'N'
  882. *
  883. IF( LOWER ) THEN
  884. *
  885. * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'L'
  886. *
  887. IF( NOTRANS ) THEN
  888. *
  889. * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
  890. * and TRANS = 'N'
  891. *
  892. CALL CTRSM( 'R', 'U', 'C', DIAG, M, K, ALPHA,
  893. $ A( 0 ), N+1, B( 0, K ), LDB )
  894. CALL CGEMM( 'N', 'N', M, K, K, -CONE, B( 0, K ),
  895. $ LDB, A( K+1 ), N+1, ALPHA, B( 0, 0 ),
  896. $ LDB )
  897. CALL CTRSM( 'R', 'L', 'N', DIAG, M, K, CONE,
  898. $ A( 1 ), N+1, B( 0, 0 ), LDB )
  899. *
  900. ELSE
  901. *
  902. * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'L',
  903. * and TRANS = 'C'
  904. *
  905. CALL CTRSM( 'R', 'L', 'C', DIAG, M, K, ALPHA,
  906. $ A( 1 ), N+1, B( 0, 0 ), LDB )
  907. CALL CGEMM( 'N', 'C', M, K, K, -CONE, B( 0, 0 ),
  908. $ LDB, A( K+1 ), N+1, ALPHA, B( 0, K ),
  909. $ LDB )
  910. CALL CTRSM( 'R', 'U', 'N', DIAG, M, K, CONE,
  911. $ A( 0 ), N+1, B( 0, K ), LDB )
  912. *
  913. END IF
  914. *
  915. ELSE
  916. *
  917. * SIDE ='R', N is even, TRANSR = 'N', and UPLO = 'U'
  918. *
  919. IF( NOTRANS ) THEN
  920. *
  921. * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
  922. * and TRANS = 'N'
  923. *
  924. CALL CTRSM( 'R', 'L', 'C', DIAG, M, K, ALPHA,
  925. $ A( K+1 ), N+1, B( 0, 0 ), LDB )
  926. CALL CGEMM( 'N', 'N', M, K, K, -CONE, B( 0, 0 ),
  927. $ LDB, A( 0 ), N+1, ALPHA, B( 0, K ),
  928. $ LDB )
  929. CALL CTRSM( 'R', 'U', 'N', DIAG, M, K, CONE,
  930. $ A( K ), N+1, B( 0, K ), LDB )
  931. *
  932. ELSE
  933. *
  934. * SIDE ='R', N is even, TRANSR = 'N', UPLO = 'U',
  935. * and TRANS = 'C'
  936. *
  937. CALL CTRSM( 'R', 'U', 'C', DIAG, M, K, ALPHA,
  938. $ A( K ), N+1, B( 0, K ), LDB )
  939. CALL CGEMM( 'N', 'C', M, K, K, -CONE, B( 0, K ),
  940. $ LDB, A( 0 ), N+1, ALPHA, B( 0, 0 ),
  941. $ LDB )
  942. CALL CTRSM( 'R', 'L', 'N', DIAG, M, K, CONE,
  943. $ A( K+1 ), N+1, B( 0, 0 ), LDB )
  944. *
  945. END IF
  946. *
  947. END IF
  948. *
  949. ELSE
  950. *
  951. * SIDE = 'R', N is even, and TRANSR = 'C'
  952. *
  953. IF( LOWER ) THEN
  954. *
  955. * SIDE ='R', N is even, TRANSR = 'C', and UPLO = 'L'
  956. *
  957. IF( NOTRANS ) THEN
  958. *
  959. * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'L',
  960. * and TRANS = 'N'
  961. *
  962. CALL CTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
  963. $ A( 0 ), K, B( 0, K ), LDB )
  964. CALL CGEMM( 'N', 'C', M, K, K, -CONE, B( 0, K ),
  965. $ LDB, A( ( K+1 )*K ), K, ALPHA,
  966. $ B( 0, 0 ), LDB )
  967. CALL CTRSM( 'R', 'U', 'C', DIAG, M, K, CONE,
  968. $ A( K ), K, B( 0, 0 ), LDB )
  969. *
  970. ELSE
  971. *
  972. * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'L',
  973. * and TRANS = 'C'
  974. *
  975. CALL CTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
  976. $ A( K ), K, B( 0, 0 ), LDB )
  977. CALL CGEMM( 'N', 'N', M, K, K, -CONE, B( 0, 0 ),
  978. $ LDB, A( ( K+1 )*K ), K, ALPHA,
  979. $ B( 0, K ), LDB )
  980. CALL CTRSM( 'R', 'L', 'C', DIAG, M, K, CONE,
  981. $ A( 0 ), K, B( 0, K ), LDB )
  982. *
  983. END IF
  984. *
  985. ELSE
  986. *
  987. * SIDE ='R', N is even, TRANSR = 'C', and UPLO = 'U'
  988. *
  989. IF( NOTRANS ) THEN
  990. *
  991. * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'U',
  992. * and TRANS = 'N'
  993. *
  994. CALL CTRSM( 'R', 'U', 'N', DIAG, M, K, ALPHA,
  995. $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
  996. CALL CGEMM( 'N', 'C', M, K, K, -CONE, B( 0, 0 ),
  997. $ LDB, A( 0 ), K, ALPHA, B( 0, K ), LDB )
  998. CALL CTRSM( 'R', 'L', 'C', DIAG, M, K, CONE,
  999. $ A( K*K ), K, B( 0, K ), LDB )
  1000. *
  1001. ELSE
  1002. *
  1003. * SIDE ='R', N is even, TRANSR = 'C', UPLO = 'U',
  1004. * and TRANS = 'C'
  1005. *
  1006. CALL CTRSM( 'R', 'L', 'N', DIAG, M, K, ALPHA,
  1007. $ A( K*K ), K, B( 0, K ), LDB )
  1008. CALL CGEMM( 'N', 'N', M, K, K, -CONE, B( 0, K ),
  1009. $ LDB, A( 0 ), K, ALPHA, B( 0, 0 ), LDB )
  1010. CALL CTRSM( 'R', 'U', 'C', DIAG, M, K, CONE,
  1011. $ A( ( K+1 )*K ), K, B( 0, 0 ), LDB )
  1012. *
  1013. END IF
  1014. *
  1015. END IF
  1016. *
  1017. END IF
  1018. *
  1019. END IF
  1020. END IF
  1021. *
  1022. RETURN
  1023. *
  1024. * End of CTFSM
  1025. *
  1026. END