ctrsm.c 19 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698
  1. /* ctrsm.f -- translated by f2c (version 20061008).
  2. You must link the resulting object file with libf2c:
  3. on Microsoft Windows system, link with libf2c.lib;
  4. on Linux or Unix systems, link with .../path/to/libf2c.a -lm
  5. or, if you install libf2c.a in a standard place, with -lf2c -lm
  6. -- in that order, at the end of the command line, as in
  7. cc *.o -lf2c -lm
  8. Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
  9. http://www.netlib.org/f2c/libf2c.zip
  10. */
  11. #include "f2c.h"
  12. #include "blaswrap.h"
  13. /* Table of constant values */
  14. static complex c_b1 = {1.f,0.f};
  15. /* Subroutine */ int ctrsm_(char *side, char *uplo, char *transa, char *diag,
  16. integer *m, integer *n, complex *alpha, complex *a, integer *lda,
  17. complex *b, integer *ldb)
  18. {
  19. /* System generated locals */
  20. integer a_dim1, a_offset, b_dim1, b_offset, i__1, i__2, i__3, i__4, i__5,
  21. i__6, i__7;
  22. complex q__1, q__2, q__3;
  23. /* Builtin functions */
  24. void c_div(complex *, complex *, complex *), r_cnjg(complex *, complex *);
  25. /* Local variables */
  26. integer i__, j, k, info;
  27. complex temp;
  28. extern logical lsame_(char *, char *);
  29. logical lside;
  30. integer nrowa;
  31. logical upper;
  32. extern /* Subroutine */ int xerbla_(char *, integer *);
  33. logical noconj, nounit;
  34. /* .. Scalar Arguments .. */
  35. /* .. */
  36. /* .. Array Arguments .. */
  37. /* .. */
  38. /* Purpose */
  39. /* ======= */
  40. /* CTRSM solves one of the matrix equations */
  41. /* op( A )*X = alpha*B, or X*op( A ) = alpha*B, */
  42. /* where alpha is a scalar, X and B are m by n matrices, A is a unit, or */
  43. /* non-unit, upper or lower triangular matrix and op( A ) is one of */
  44. /* op( A ) = A or op( A ) = A' or op( A ) = conjg( A' ). */
  45. /* The matrix X is overwritten on B. */
  46. /* Arguments */
  47. /* ========== */
  48. /* SIDE - CHARACTER*1. */
  49. /* On entry, SIDE specifies whether op( A ) appears on the left */
  50. /* or right of X as follows: */
  51. /* SIDE = 'L' or 'l' op( A )*X = alpha*B. */
  52. /* SIDE = 'R' or 'r' X*op( A ) = alpha*B. */
  53. /* Unchanged on exit. */
  54. /* UPLO - CHARACTER*1. */
  55. /* On entry, UPLO specifies whether the matrix A is an upper or */
  56. /* lower triangular matrix as follows: */
  57. /* UPLO = 'U' or 'u' A is an upper triangular matrix. */
  58. /* UPLO = 'L' or 'l' A is a lower triangular matrix. */
  59. /* Unchanged on exit. */
  60. /* TRANSA - CHARACTER*1. */
  61. /* On entry, TRANSA specifies the form of op( A ) to be used in */
  62. /* the matrix multiplication as follows: */
  63. /* TRANSA = 'N' or 'n' op( A ) = A. */
  64. /* TRANSA = 'T' or 't' op( A ) = A'. */
  65. /* TRANSA = 'C' or 'c' op( A ) = conjg( A' ). */
  66. /* Unchanged on exit. */
  67. /* DIAG - CHARACTER*1. */
  68. /* On entry, DIAG specifies whether or not A is unit triangular */
  69. /* as follows: */
  70. /* DIAG = 'U' or 'u' A is assumed to be unit triangular. */
  71. /* DIAG = 'N' or 'n' A is not assumed to be unit */
  72. /* triangular. */
  73. /* Unchanged on exit. */
  74. /* M - INTEGER. */
  75. /* On entry, M specifies the number of rows of B. M must be at */
  76. /* least zero. */
  77. /* Unchanged on exit. */
  78. /* N - INTEGER. */
  79. /* On entry, N specifies the number of columns of B. N must be */
  80. /* at least zero. */
  81. /* Unchanged on exit. */
  82. /* ALPHA - COMPLEX . */
  83. /* On entry, ALPHA specifies the scalar alpha. When alpha is */
  84. /* zero then A is not referenced and B need not be set before */
  85. /* entry. */
  86. /* Unchanged on exit. */
  87. /* A - COMPLEX array of DIMENSION ( LDA, k ), where k is m */
  88. /* when SIDE = 'L' or 'l' and is n when SIDE = 'R' or 'r'. */
  89. /* Before entry with UPLO = 'U' or 'u', the leading k by k */
  90. /* upper triangular part of the array A must contain the upper */
  91. /* triangular matrix and the strictly lower triangular part of */
  92. /* A is not referenced. */
  93. /* Before entry with UPLO = 'L' or 'l', the leading k by k */
  94. /* lower triangular part of the array A must contain the lower */
  95. /* triangular matrix and the strictly upper triangular part of */
  96. /* A is not referenced. */
  97. /* Note that when DIAG = 'U' or 'u', the diagonal elements of */
  98. /* A are not referenced either, but are assumed to be unity. */
  99. /* Unchanged on exit. */
  100. /* LDA - INTEGER. */
  101. /* On entry, LDA specifies the first dimension of A as declared */
  102. /* in the calling (sub) program. When SIDE = 'L' or 'l' then */
  103. /* LDA must be at least max( 1, m ), when SIDE = 'R' or 'r' */
  104. /* then LDA must be at least max( 1, n ). */
  105. /* Unchanged on exit. */
  106. /* B - COMPLEX array of DIMENSION ( LDB, n ). */
  107. /* Before entry, the leading m by n part of the array B must */
  108. /* contain the right-hand side matrix B, and on exit is */
  109. /* overwritten by the solution matrix X. */
  110. /* LDB - INTEGER. */
  111. /* On entry, LDB specifies the first dimension of B as declared */
  112. /* in the calling (sub) program. LDB must be at least */
  113. /* max( 1, m ). */
  114. /* Unchanged on exit. */
  115. /* Level 3 Blas routine. */
  116. /* -- Written on 8-February-1989. */
  117. /* Jack Dongarra, Argonne National Laboratory. */
  118. /* Iain Duff, AERE Harwell. */
  119. /* Jeremy Du Croz, Numerical Algorithms Group Ltd. */
  120. /* Sven Hammarling, Numerical Algorithms Group Ltd. */
  121. /* .. External Functions .. */
  122. /* .. */
  123. /* .. External Subroutines .. */
  124. /* .. */
  125. /* .. Intrinsic Functions .. */
  126. /* .. */
  127. /* .. Local Scalars .. */
  128. /* .. */
  129. /* .. Parameters .. */
  130. /* .. */
  131. /* Test the input parameters. */
  132. /* Parameter adjustments */
  133. a_dim1 = *lda;
  134. a_offset = 1 + a_dim1;
  135. a -= a_offset;
  136. b_dim1 = *ldb;
  137. b_offset = 1 + b_dim1;
  138. b -= b_offset;
  139. /* Function Body */
  140. lside = lsame_(side, "L");
  141. if (lside) {
  142. nrowa = *m;
  143. } else {
  144. nrowa = *n;
  145. }
  146. noconj = lsame_(transa, "T");
  147. nounit = lsame_(diag, "N");
  148. upper = lsame_(uplo, "U");
  149. info = 0;
  150. if (! lside && ! lsame_(side, "R")) {
  151. info = 1;
  152. } else if (! upper && ! lsame_(uplo, "L")) {
  153. info = 2;
  154. } else if (! lsame_(transa, "N") && ! lsame_(transa,
  155. "T") && ! lsame_(transa, "C")) {
  156. info = 3;
  157. } else if (! lsame_(diag, "U") && ! lsame_(diag,
  158. "N")) {
  159. info = 4;
  160. } else if (*m < 0) {
  161. info = 5;
  162. } else if (*n < 0) {
  163. info = 6;
  164. } else if (*lda < max(1,nrowa)) {
  165. info = 9;
  166. } else if (*ldb < max(1,*m)) {
  167. info = 11;
  168. }
  169. if (info != 0) {
  170. xerbla_("CTRSM ", &info);
  171. return 0;
  172. }
  173. /* Quick return if possible. */
  174. if (*m == 0 || *n == 0) {
  175. return 0;
  176. }
  177. /* And when alpha.eq.zero. */
  178. if (alpha->r == 0.f && alpha->i == 0.f) {
  179. i__1 = *n;
  180. for (j = 1; j <= i__1; ++j) {
  181. i__2 = *m;
  182. for (i__ = 1; i__ <= i__2; ++i__) {
  183. i__3 = i__ + j * b_dim1;
  184. b[i__3].r = 0.f, b[i__3].i = 0.f;
  185. /* L10: */
  186. }
  187. /* L20: */
  188. }
  189. return 0;
  190. }
  191. /* Start the operations. */
  192. if (lside) {
  193. if (lsame_(transa, "N")) {
  194. /* Form B := alpha*inv( A )*B. */
  195. if (upper) {
  196. i__1 = *n;
  197. for (j = 1; j <= i__1; ++j) {
  198. if (alpha->r != 1.f || alpha->i != 0.f) {
  199. i__2 = *m;
  200. for (i__ = 1; i__ <= i__2; ++i__) {
  201. i__3 = i__ + j * b_dim1;
  202. i__4 = i__ + j * b_dim1;
  203. q__1.r = alpha->r * b[i__4].r - alpha->i * b[i__4]
  204. .i, q__1.i = alpha->r * b[i__4].i +
  205. alpha->i * b[i__4].r;
  206. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  207. /* L30: */
  208. }
  209. }
  210. for (k = *m; k >= 1; --k) {
  211. i__2 = k + j * b_dim1;
  212. if (b[i__2].r != 0.f || b[i__2].i != 0.f) {
  213. if (nounit) {
  214. i__2 = k + j * b_dim1;
  215. c_div(&q__1, &b[k + j * b_dim1], &a[k + k *
  216. a_dim1]);
  217. b[i__2].r = q__1.r, b[i__2].i = q__1.i;
  218. }
  219. i__2 = k - 1;
  220. for (i__ = 1; i__ <= i__2; ++i__) {
  221. i__3 = i__ + j * b_dim1;
  222. i__4 = i__ + j * b_dim1;
  223. i__5 = k + j * b_dim1;
  224. i__6 = i__ + k * a_dim1;
  225. q__2.r = b[i__5].r * a[i__6].r - b[i__5].i *
  226. a[i__6].i, q__2.i = b[i__5].r * a[
  227. i__6].i + b[i__5].i * a[i__6].r;
  228. q__1.r = b[i__4].r - q__2.r, q__1.i = b[i__4]
  229. .i - q__2.i;
  230. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  231. /* L40: */
  232. }
  233. }
  234. /* L50: */
  235. }
  236. /* L60: */
  237. }
  238. } else {
  239. i__1 = *n;
  240. for (j = 1; j <= i__1; ++j) {
  241. if (alpha->r != 1.f || alpha->i != 0.f) {
  242. i__2 = *m;
  243. for (i__ = 1; i__ <= i__2; ++i__) {
  244. i__3 = i__ + j * b_dim1;
  245. i__4 = i__ + j * b_dim1;
  246. q__1.r = alpha->r * b[i__4].r - alpha->i * b[i__4]
  247. .i, q__1.i = alpha->r * b[i__4].i +
  248. alpha->i * b[i__4].r;
  249. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  250. /* L70: */
  251. }
  252. }
  253. i__2 = *m;
  254. for (k = 1; k <= i__2; ++k) {
  255. i__3 = k + j * b_dim1;
  256. if (b[i__3].r != 0.f || b[i__3].i != 0.f) {
  257. if (nounit) {
  258. i__3 = k + j * b_dim1;
  259. c_div(&q__1, &b[k + j * b_dim1], &a[k + k *
  260. a_dim1]);
  261. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  262. }
  263. i__3 = *m;
  264. for (i__ = k + 1; i__ <= i__3; ++i__) {
  265. i__4 = i__ + j * b_dim1;
  266. i__5 = i__ + j * b_dim1;
  267. i__6 = k + j * b_dim1;
  268. i__7 = i__ + k * a_dim1;
  269. q__2.r = b[i__6].r * a[i__7].r - b[i__6].i *
  270. a[i__7].i, q__2.i = b[i__6].r * a[
  271. i__7].i + b[i__6].i * a[i__7].r;
  272. q__1.r = b[i__5].r - q__2.r, q__1.i = b[i__5]
  273. .i - q__2.i;
  274. b[i__4].r = q__1.r, b[i__4].i = q__1.i;
  275. /* L80: */
  276. }
  277. }
  278. /* L90: */
  279. }
  280. /* L100: */
  281. }
  282. }
  283. } else {
  284. /* Form B := alpha*inv( A' )*B */
  285. /* or B := alpha*inv( conjg( A' ) )*B. */
  286. if (upper) {
  287. i__1 = *n;
  288. for (j = 1; j <= i__1; ++j) {
  289. i__2 = *m;
  290. for (i__ = 1; i__ <= i__2; ++i__) {
  291. i__3 = i__ + j * b_dim1;
  292. q__1.r = alpha->r * b[i__3].r - alpha->i * b[i__3].i,
  293. q__1.i = alpha->r * b[i__3].i + alpha->i * b[
  294. i__3].r;
  295. temp.r = q__1.r, temp.i = q__1.i;
  296. if (noconj) {
  297. i__3 = i__ - 1;
  298. for (k = 1; k <= i__3; ++k) {
  299. i__4 = k + i__ * a_dim1;
  300. i__5 = k + j * b_dim1;
  301. q__2.r = a[i__4].r * b[i__5].r - a[i__4].i *
  302. b[i__5].i, q__2.i = a[i__4].r * b[
  303. i__5].i + a[i__4].i * b[i__5].r;
  304. q__1.r = temp.r - q__2.r, q__1.i = temp.i -
  305. q__2.i;
  306. temp.r = q__1.r, temp.i = q__1.i;
  307. /* L110: */
  308. }
  309. if (nounit) {
  310. c_div(&q__1, &temp, &a[i__ + i__ * a_dim1]);
  311. temp.r = q__1.r, temp.i = q__1.i;
  312. }
  313. } else {
  314. i__3 = i__ - 1;
  315. for (k = 1; k <= i__3; ++k) {
  316. r_cnjg(&q__3, &a[k + i__ * a_dim1]);
  317. i__4 = k + j * b_dim1;
  318. q__2.r = q__3.r * b[i__4].r - q__3.i * b[i__4]
  319. .i, q__2.i = q__3.r * b[i__4].i +
  320. q__3.i * b[i__4].r;
  321. q__1.r = temp.r - q__2.r, q__1.i = temp.i -
  322. q__2.i;
  323. temp.r = q__1.r, temp.i = q__1.i;
  324. /* L120: */
  325. }
  326. if (nounit) {
  327. r_cnjg(&q__2, &a[i__ + i__ * a_dim1]);
  328. c_div(&q__1, &temp, &q__2);
  329. temp.r = q__1.r, temp.i = q__1.i;
  330. }
  331. }
  332. i__3 = i__ + j * b_dim1;
  333. b[i__3].r = temp.r, b[i__3].i = temp.i;
  334. /* L130: */
  335. }
  336. /* L140: */
  337. }
  338. } else {
  339. i__1 = *n;
  340. for (j = 1; j <= i__1; ++j) {
  341. for (i__ = *m; i__ >= 1; --i__) {
  342. i__2 = i__ + j * b_dim1;
  343. q__1.r = alpha->r * b[i__2].r - alpha->i * b[i__2].i,
  344. q__1.i = alpha->r * b[i__2].i + alpha->i * b[
  345. i__2].r;
  346. temp.r = q__1.r, temp.i = q__1.i;
  347. if (noconj) {
  348. i__2 = *m;
  349. for (k = i__ + 1; k <= i__2; ++k) {
  350. i__3 = k + i__ * a_dim1;
  351. i__4 = k + j * b_dim1;
  352. q__2.r = a[i__3].r * b[i__4].r - a[i__3].i *
  353. b[i__4].i, q__2.i = a[i__3].r * b[
  354. i__4].i + a[i__3].i * b[i__4].r;
  355. q__1.r = temp.r - q__2.r, q__1.i = temp.i -
  356. q__2.i;
  357. temp.r = q__1.r, temp.i = q__1.i;
  358. /* L150: */
  359. }
  360. if (nounit) {
  361. c_div(&q__1, &temp, &a[i__ + i__ * a_dim1]);
  362. temp.r = q__1.r, temp.i = q__1.i;
  363. }
  364. } else {
  365. i__2 = *m;
  366. for (k = i__ + 1; k <= i__2; ++k) {
  367. r_cnjg(&q__3, &a[k + i__ * a_dim1]);
  368. i__3 = k + j * b_dim1;
  369. q__2.r = q__3.r * b[i__3].r - q__3.i * b[i__3]
  370. .i, q__2.i = q__3.r * b[i__3].i +
  371. q__3.i * b[i__3].r;
  372. q__1.r = temp.r - q__2.r, q__1.i = temp.i -
  373. q__2.i;
  374. temp.r = q__1.r, temp.i = q__1.i;
  375. /* L160: */
  376. }
  377. if (nounit) {
  378. r_cnjg(&q__2, &a[i__ + i__ * a_dim1]);
  379. c_div(&q__1, &temp, &q__2);
  380. temp.r = q__1.r, temp.i = q__1.i;
  381. }
  382. }
  383. i__2 = i__ + j * b_dim1;
  384. b[i__2].r = temp.r, b[i__2].i = temp.i;
  385. /* L170: */
  386. }
  387. /* L180: */
  388. }
  389. }
  390. }
  391. } else {
  392. if (lsame_(transa, "N")) {
  393. /* Form B := alpha*B*inv( A ). */
  394. if (upper) {
  395. i__1 = *n;
  396. for (j = 1; j <= i__1; ++j) {
  397. if (alpha->r != 1.f || alpha->i != 0.f) {
  398. i__2 = *m;
  399. for (i__ = 1; i__ <= i__2; ++i__) {
  400. i__3 = i__ + j * b_dim1;
  401. i__4 = i__ + j * b_dim1;
  402. q__1.r = alpha->r * b[i__4].r - alpha->i * b[i__4]
  403. .i, q__1.i = alpha->r * b[i__4].i +
  404. alpha->i * b[i__4].r;
  405. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  406. /* L190: */
  407. }
  408. }
  409. i__2 = j - 1;
  410. for (k = 1; k <= i__2; ++k) {
  411. i__3 = k + j * a_dim1;
  412. if (a[i__3].r != 0.f || a[i__3].i != 0.f) {
  413. i__3 = *m;
  414. for (i__ = 1; i__ <= i__3; ++i__) {
  415. i__4 = i__ + j * b_dim1;
  416. i__5 = i__ + j * b_dim1;
  417. i__6 = k + j * a_dim1;
  418. i__7 = i__ + k * b_dim1;
  419. q__2.r = a[i__6].r * b[i__7].r - a[i__6].i *
  420. b[i__7].i, q__2.i = a[i__6].r * b[
  421. i__7].i + a[i__6].i * b[i__7].r;
  422. q__1.r = b[i__5].r - q__2.r, q__1.i = b[i__5]
  423. .i - q__2.i;
  424. b[i__4].r = q__1.r, b[i__4].i = q__1.i;
  425. /* L200: */
  426. }
  427. }
  428. /* L210: */
  429. }
  430. if (nounit) {
  431. c_div(&q__1, &c_b1, &a[j + j * a_dim1]);
  432. temp.r = q__1.r, temp.i = q__1.i;
  433. i__2 = *m;
  434. for (i__ = 1; i__ <= i__2; ++i__) {
  435. i__3 = i__ + j * b_dim1;
  436. i__4 = i__ + j * b_dim1;
  437. q__1.r = temp.r * b[i__4].r - temp.i * b[i__4].i,
  438. q__1.i = temp.r * b[i__4].i + temp.i * b[
  439. i__4].r;
  440. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  441. /* L220: */
  442. }
  443. }
  444. /* L230: */
  445. }
  446. } else {
  447. for (j = *n; j >= 1; --j) {
  448. if (alpha->r != 1.f || alpha->i != 0.f) {
  449. i__1 = *m;
  450. for (i__ = 1; i__ <= i__1; ++i__) {
  451. i__2 = i__ + j * b_dim1;
  452. i__3 = i__ + j * b_dim1;
  453. q__1.r = alpha->r * b[i__3].r - alpha->i * b[i__3]
  454. .i, q__1.i = alpha->r * b[i__3].i +
  455. alpha->i * b[i__3].r;
  456. b[i__2].r = q__1.r, b[i__2].i = q__1.i;
  457. /* L240: */
  458. }
  459. }
  460. i__1 = *n;
  461. for (k = j + 1; k <= i__1; ++k) {
  462. i__2 = k + j * a_dim1;
  463. if (a[i__2].r != 0.f || a[i__2].i != 0.f) {
  464. i__2 = *m;
  465. for (i__ = 1; i__ <= i__2; ++i__) {
  466. i__3 = i__ + j * b_dim1;
  467. i__4 = i__ + j * b_dim1;
  468. i__5 = k + j * a_dim1;
  469. i__6 = i__ + k * b_dim1;
  470. q__2.r = a[i__5].r * b[i__6].r - a[i__5].i *
  471. b[i__6].i, q__2.i = a[i__5].r * b[
  472. i__6].i + a[i__5].i * b[i__6].r;
  473. q__1.r = b[i__4].r - q__2.r, q__1.i = b[i__4]
  474. .i - q__2.i;
  475. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  476. /* L250: */
  477. }
  478. }
  479. /* L260: */
  480. }
  481. if (nounit) {
  482. c_div(&q__1, &c_b1, &a[j + j * a_dim1]);
  483. temp.r = q__1.r, temp.i = q__1.i;
  484. i__1 = *m;
  485. for (i__ = 1; i__ <= i__1; ++i__) {
  486. i__2 = i__ + j * b_dim1;
  487. i__3 = i__ + j * b_dim1;
  488. q__1.r = temp.r * b[i__3].r - temp.i * b[i__3].i,
  489. q__1.i = temp.r * b[i__3].i + temp.i * b[
  490. i__3].r;
  491. b[i__2].r = q__1.r, b[i__2].i = q__1.i;
  492. /* L270: */
  493. }
  494. }
  495. /* L280: */
  496. }
  497. }
  498. } else {
  499. /* Form B := alpha*B*inv( A' ) */
  500. /* or B := alpha*B*inv( conjg( A' ) ). */
  501. if (upper) {
  502. for (k = *n; k >= 1; --k) {
  503. if (nounit) {
  504. if (noconj) {
  505. c_div(&q__1, &c_b1, &a[k + k * a_dim1]);
  506. temp.r = q__1.r, temp.i = q__1.i;
  507. } else {
  508. r_cnjg(&q__2, &a[k + k * a_dim1]);
  509. c_div(&q__1, &c_b1, &q__2);
  510. temp.r = q__1.r, temp.i = q__1.i;
  511. }
  512. i__1 = *m;
  513. for (i__ = 1; i__ <= i__1; ++i__) {
  514. i__2 = i__ + k * b_dim1;
  515. i__3 = i__ + k * b_dim1;
  516. q__1.r = temp.r * b[i__3].r - temp.i * b[i__3].i,
  517. q__1.i = temp.r * b[i__3].i + temp.i * b[
  518. i__3].r;
  519. b[i__2].r = q__1.r, b[i__2].i = q__1.i;
  520. /* L290: */
  521. }
  522. }
  523. i__1 = k - 1;
  524. for (j = 1; j <= i__1; ++j) {
  525. i__2 = j + k * a_dim1;
  526. if (a[i__2].r != 0.f || a[i__2].i != 0.f) {
  527. if (noconj) {
  528. i__2 = j + k * a_dim1;
  529. temp.r = a[i__2].r, temp.i = a[i__2].i;
  530. } else {
  531. r_cnjg(&q__1, &a[j + k * a_dim1]);
  532. temp.r = q__1.r, temp.i = q__1.i;
  533. }
  534. i__2 = *m;
  535. for (i__ = 1; i__ <= i__2; ++i__) {
  536. i__3 = i__ + j * b_dim1;
  537. i__4 = i__ + j * b_dim1;
  538. i__5 = i__ + k * b_dim1;
  539. q__2.r = temp.r * b[i__5].r - temp.i * b[i__5]
  540. .i, q__2.i = temp.r * b[i__5].i +
  541. temp.i * b[i__5].r;
  542. q__1.r = b[i__4].r - q__2.r, q__1.i = b[i__4]
  543. .i - q__2.i;
  544. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  545. /* L300: */
  546. }
  547. }
  548. /* L310: */
  549. }
  550. if (alpha->r != 1.f || alpha->i != 0.f) {
  551. i__1 = *m;
  552. for (i__ = 1; i__ <= i__1; ++i__) {
  553. i__2 = i__ + k * b_dim1;
  554. i__3 = i__ + k * b_dim1;
  555. q__1.r = alpha->r * b[i__3].r - alpha->i * b[i__3]
  556. .i, q__1.i = alpha->r * b[i__3].i +
  557. alpha->i * b[i__3].r;
  558. b[i__2].r = q__1.r, b[i__2].i = q__1.i;
  559. /* L320: */
  560. }
  561. }
  562. /* L330: */
  563. }
  564. } else {
  565. i__1 = *n;
  566. for (k = 1; k <= i__1; ++k) {
  567. if (nounit) {
  568. if (noconj) {
  569. c_div(&q__1, &c_b1, &a[k + k * a_dim1]);
  570. temp.r = q__1.r, temp.i = q__1.i;
  571. } else {
  572. r_cnjg(&q__2, &a[k + k * a_dim1]);
  573. c_div(&q__1, &c_b1, &q__2);
  574. temp.r = q__1.r, temp.i = q__1.i;
  575. }
  576. i__2 = *m;
  577. for (i__ = 1; i__ <= i__2; ++i__) {
  578. i__3 = i__ + k * b_dim1;
  579. i__4 = i__ + k * b_dim1;
  580. q__1.r = temp.r * b[i__4].r - temp.i * b[i__4].i,
  581. q__1.i = temp.r * b[i__4].i + temp.i * b[
  582. i__4].r;
  583. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  584. /* L340: */
  585. }
  586. }
  587. i__2 = *n;
  588. for (j = k + 1; j <= i__2; ++j) {
  589. i__3 = j + k * a_dim1;
  590. if (a[i__3].r != 0.f || a[i__3].i != 0.f) {
  591. if (noconj) {
  592. i__3 = j + k * a_dim1;
  593. temp.r = a[i__3].r, temp.i = a[i__3].i;
  594. } else {
  595. r_cnjg(&q__1, &a[j + k * a_dim1]);
  596. temp.r = q__1.r, temp.i = q__1.i;
  597. }
  598. i__3 = *m;
  599. for (i__ = 1; i__ <= i__3; ++i__) {
  600. i__4 = i__ + j * b_dim1;
  601. i__5 = i__ + j * b_dim1;
  602. i__6 = i__ + k * b_dim1;
  603. q__2.r = temp.r * b[i__6].r - temp.i * b[i__6]
  604. .i, q__2.i = temp.r * b[i__6].i +
  605. temp.i * b[i__6].r;
  606. q__1.r = b[i__5].r - q__2.r, q__1.i = b[i__5]
  607. .i - q__2.i;
  608. b[i__4].r = q__1.r, b[i__4].i = q__1.i;
  609. /* L350: */
  610. }
  611. }
  612. /* L360: */
  613. }
  614. if (alpha->r != 1.f || alpha->i != 0.f) {
  615. i__2 = *m;
  616. for (i__ = 1; i__ <= i__2; ++i__) {
  617. i__3 = i__ + k * b_dim1;
  618. i__4 = i__ + k * b_dim1;
  619. q__1.r = alpha->r * b[i__4].r - alpha->i * b[i__4]
  620. .i, q__1.i = alpha->r * b[i__4].i +
  621. alpha->i * b[i__4].r;
  622. b[i__3].r = q__1.r, b[i__3].i = q__1.i;
  623. /* L370: */
  624. }
  625. }
  626. /* L380: */
  627. }
  628. }
  629. }
  630. }
  631. return 0;
  632. /* End of CTRSM . */
  633. } /* ctrsm_ */