00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015
00016
00017
00018 static complex c_b1 = {1.f,0.f};
00019
00020 int ctrsm_(char *side, char *uplo, char *transa, char *diag,
00021 integer *m, integer *n, complex *alpha, complex *a, integer *lda,
00022 complex *b, integer *ldb)
00023 {
00024
00025 integer a_dim1, a_offset, b_dim1, b_offset, i__1, i__2, i__3, i__4, i__5,
00026 i__6, i__7;
00027 complex q__1, q__2, q__3;
00028
00029
00030 void c_div(complex *, complex *, complex *), r_cnjg(complex *, complex *);
00031
00032
00033 integer i__, j, k, info;
00034 complex temp;
00035 extern logical lsame_(char *, char *);
00036 logical lside;
00037 integer nrowa;
00038 logical upper;
00039 extern int xerbla_(char *, integer *);
00040 logical noconj, nounit;
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179 a_dim1 = *lda;
00180 a_offset = 1 + a_dim1;
00181 a -= a_offset;
00182 b_dim1 = *ldb;
00183 b_offset = 1 + b_dim1;
00184 b -= b_offset;
00185
00186
00187 lside = lsame_(side, "L");
00188 if (lside) {
00189 nrowa = *m;
00190 } else {
00191 nrowa = *n;
00192 }
00193 noconj = lsame_(transa, "T");
00194 nounit = lsame_(diag, "N");
00195 upper = lsame_(uplo, "U");
00196
00197 info = 0;
00198 if (! lside && ! lsame_(side, "R")) {
00199 info = 1;
00200 } else if (! upper && ! lsame_(uplo, "L")) {
00201 info = 2;
00202 } else if (! lsame_(transa, "N") && ! lsame_(transa,
00203 "T") && ! lsame_(transa, "C")) {
00204 info = 3;
00205 } else if (! lsame_(diag, "U") && ! lsame_(diag,
00206 "N")) {
00207 info = 4;
00208 } else if (*m < 0) {
00209 info = 5;
00210 } else if (*n < 0) {
00211 info = 6;
00212 } else if (*lda < max(1,nrowa)) {
00213 info = 9;
00214 } else if (*ldb < max(1,*m)) {
00215 info = 11;
00216 }
00217 if (info != 0) {
00218 xerbla_("CTRSM ", &info);
00219 return 0;
00220 }
00221
00222
00223
00224 if (*m == 0 || *n == 0) {
00225 return 0;
00226 }
00227
00228
00229
00230 if (alpha->r == 0.f && alpha->i == 0.f) {
00231 i__1 = *n;
00232 for (j = 1; j <= i__1; ++j) {
00233 i__2 = *m;
00234 for (i__ = 1; i__ <= i__2; ++i__) {
00235 i__3 = i__ + j * b_dim1;
00236 b[i__3].r = 0.f, b[i__3].i = 0.f;
00237
00238 }
00239
00240 }
00241 return 0;
00242 }
00243
00244
00245
00246 if (lside) {
00247 if (lsame_(transa, "N")) {
00248
00249
00250
00251 if (upper) {
00252 i__1 = *n;
00253 for (j = 1; j <= i__1; ++j) {
00254 if (alpha->r != 1.f || alpha->i != 0.f) {
00255 i__2 = *m;
00256 for (i__ = 1; i__ <= i__2; ++i__) {
00257 i__3 = i__ + j * b_dim1;
00258 i__4 = i__ + j * b_dim1;
00259 q__1.r = alpha->r * b[i__4].r - alpha->i * b[i__4]
00260 .i, q__1.i = alpha->r * b[i__4].i +
00261 alpha->i * b[i__4].r;
00262 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00263
00264 }
00265 }
00266 for (k = *m; k >= 1; --k) {
00267 i__2 = k + j * b_dim1;
00268 if (b[i__2].r != 0.f || b[i__2].i != 0.f) {
00269 if (nounit) {
00270 i__2 = k + j * b_dim1;
00271 c_div(&q__1, &b[k + j * b_dim1], &a[k + k *
00272 a_dim1]);
00273 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00274 }
00275 i__2 = k - 1;
00276 for (i__ = 1; i__ <= i__2; ++i__) {
00277 i__3 = i__ + j * b_dim1;
00278 i__4 = i__ + j * b_dim1;
00279 i__5 = k + j * b_dim1;
00280 i__6 = i__ + k * a_dim1;
00281 q__2.r = b[i__5].r * a[i__6].r - b[i__5].i *
00282 a[i__6].i, q__2.i = b[i__5].r * a[
00283 i__6].i + b[i__5].i * a[i__6].r;
00284 q__1.r = b[i__4].r - q__2.r, q__1.i = b[i__4]
00285 .i - q__2.i;
00286 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00287
00288 }
00289 }
00290
00291 }
00292
00293 }
00294 } else {
00295 i__1 = *n;
00296 for (j = 1; j <= i__1; ++j) {
00297 if (alpha->r != 1.f || alpha->i != 0.f) {
00298 i__2 = *m;
00299 for (i__ = 1; i__ <= i__2; ++i__) {
00300 i__3 = i__ + j * b_dim1;
00301 i__4 = i__ + j * b_dim1;
00302 q__1.r = alpha->r * b[i__4].r - alpha->i * b[i__4]
00303 .i, q__1.i = alpha->r * b[i__4].i +
00304 alpha->i * b[i__4].r;
00305 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00306
00307 }
00308 }
00309 i__2 = *m;
00310 for (k = 1; k <= i__2; ++k) {
00311 i__3 = k + j * b_dim1;
00312 if (b[i__3].r != 0.f || b[i__3].i != 0.f) {
00313 if (nounit) {
00314 i__3 = k + j * b_dim1;
00315 c_div(&q__1, &b[k + j * b_dim1], &a[k + k *
00316 a_dim1]);
00317 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00318 }
00319 i__3 = *m;
00320 for (i__ = k + 1; i__ <= i__3; ++i__) {
00321 i__4 = i__ + j * b_dim1;
00322 i__5 = i__ + j * b_dim1;
00323 i__6 = k + j * b_dim1;
00324 i__7 = i__ + k * a_dim1;
00325 q__2.r = b[i__6].r * a[i__7].r - b[i__6].i *
00326 a[i__7].i, q__2.i = b[i__6].r * a[
00327 i__7].i + b[i__6].i * a[i__7].r;
00328 q__1.r = b[i__5].r - q__2.r, q__1.i = b[i__5]
00329 .i - q__2.i;
00330 b[i__4].r = q__1.r, b[i__4].i = q__1.i;
00331
00332 }
00333 }
00334
00335 }
00336
00337 }
00338 }
00339 } else {
00340
00341
00342
00343
00344 if (upper) {
00345 i__1 = *n;
00346 for (j = 1; j <= i__1; ++j) {
00347 i__2 = *m;
00348 for (i__ = 1; i__ <= i__2; ++i__) {
00349 i__3 = i__ + j * b_dim1;
00350 q__1.r = alpha->r * b[i__3].r - alpha->i * b[i__3].i,
00351 q__1.i = alpha->r * b[i__3].i + alpha->i * b[
00352 i__3].r;
00353 temp.r = q__1.r, temp.i = q__1.i;
00354 if (noconj) {
00355 i__3 = i__ - 1;
00356 for (k = 1; k <= i__3; ++k) {
00357 i__4 = k + i__ * a_dim1;
00358 i__5 = k + j * b_dim1;
00359 q__2.r = a[i__4].r * b[i__5].r - a[i__4].i *
00360 b[i__5].i, q__2.i = a[i__4].r * b[
00361 i__5].i + a[i__4].i * b[i__5].r;
00362 q__1.r = temp.r - q__2.r, q__1.i = temp.i -
00363 q__2.i;
00364 temp.r = q__1.r, temp.i = q__1.i;
00365
00366 }
00367 if (nounit) {
00368 c_div(&q__1, &temp, &a[i__ + i__ * a_dim1]);
00369 temp.r = q__1.r, temp.i = q__1.i;
00370 }
00371 } else {
00372 i__3 = i__ - 1;
00373 for (k = 1; k <= i__3; ++k) {
00374 r_cnjg(&q__3, &a[k + i__ * a_dim1]);
00375 i__4 = k + j * b_dim1;
00376 q__2.r = q__3.r * b[i__4].r - q__3.i * b[i__4]
00377 .i, q__2.i = q__3.r * b[i__4].i +
00378 q__3.i * b[i__4].r;
00379 q__1.r = temp.r - q__2.r, q__1.i = temp.i -
00380 q__2.i;
00381 temp.r = q__1.r, temp.i = q__1.i;
00382
00383 }
00384 if (nounit) {
00385 r_cnjg(&q__2, &a[i__ + i__ * a_dim1]);
00386 c_div(&q__1, &temp, &q__2);
00387 temp.r = q__1.r, temp.i = q__1.i;
00388 }
00389 }
00390 i__3 = i__ + j * b_dim1;
00391 b[i__3].r = temp.r, b[i__3].i = temp.i;
00392
00393 }
00394
00395 }
00396 } else {
00397 i__1 = *n;
00398 for (j = 1; j <= i__1; ++j) {
00399 for (i__ = *m; i__ >= 1; --i__) {
00400 i__2 = i__ + j * b_dim1;
00401 q__1.r = alpha->r * b[i__2].r - alpha->i * b[i__2].i,
00402 q__1.i = alpha->r * b[i__2].i + alpha->i * b[
00403 i__2].r;
00404 temp.r = q__1.r, temp.i = q__1.i;
00405 if (noconj) {
00406 i__2 = *m;
00407 for (k = i__ + 1; k <= i__2; ++k) {
00408 i__3 = k + i__ * a_dim1;
00409 i__4 = k + j * b_dim1;
00410 q__2.r = a[i__3].r * b[i__4].r - a[i__3].i *
00411 b[i__4].i, q__2.i = a[i__3].r * b[
00412 i__4].i + a[i__3].i * b[i__4].r;
00413 q__1.r = temp.r - q__2.r, q__1.i = temp.i -
00414 q__2.i;
00415 temp.r = q__1.r, temp.i = q__1.i;
00416
00417 }
00418 if (nounit) {
00419 c_div(&q__1, &temp, &a[i__ + i__ * a_dim1]);
00420 temp.r = q__1.r, temp.i = q__1.i;
00421 }
00422 } else {
00423 i__2 = *m;
00424 for (k = i__ + 1; k <= i__2; ++k) {
00425 r_cnjg(&q__3, &a[k + i__ * a_dim1]);
00426 i__3 = k + j * b_dim1;
00427 q__2.r = q__3.r * b[i__3].r - q__3.i * b[i__3]
00428 .i, q__2.i = q__3.r * b[i__3].i +
00429 q__3.i * b[i__3].r;
00430 q__1.r = temp.r - q__2.r, q__1.i = temp.i -
00431 q__2.i;
00432 temp.r = q__1.r, temp.i = q__1.i;
00433
00434 }
00435 if (nounit) {
00436 r_cnjg(&q__2, &a[i__ + i__ * a_dim1]);
00437 c_div(&q__1, &temp, &q__2);
00438 temp.r = q__1.r, temp.i = q__1.i;
00439 }
00440 }
00441 i__2 = i__ + j * b_dim1;
00442 b[i__2].r = temp.r, b[i__2].i = temp.i;
00443
00444 }
00445
00446 }
00447 }
00448 }
00449 } else {
00450 if (lsame_(transa, "N")) {
00451
00452
00453
00454 if (upper) {
00455 i__1 = *n;
00456 for (j = 1; j <= i__1; ++j) {
00457 if (alpha->r != 1.f || alpha->i != 0.f) {
00458 i__2 = *m;
00459 for (i__ = 1; i__ <= i__2; ++i__) {
00460 i__3 = i__ + j * b_dim1;
00461 i__4 = i__ + j * b_dim1;
00462 q__1.r = alpha->r * b[i__4].r - alpha->i * b[i__4]
00463 .i, q__1.i = alpha->r * b[i__4].i +
00464 alpha->i * b[i__4].r;
00465 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00466
00467 }
00468 }
00469 i__2 = j - 1;
00470 for (k = 1; k <= i__2; ++k) {
00471 i__3 = k + j * a_dim1;
00472 if (a[i__3].r != 0.f || a[i__3].i != 0.f) {
00473 i__3 = *m;
00474 for (i__ = 1; i__ <= i__3; ++i__) {
00475 i__4 = i__ + j * b_dim1;
00476 i__5 = i__ + j * b_dim1;
00477 i__6 = k + j * a_dim1;
00478 i__7 = i__ + k * b_dim1;
00479 q__2.r = a[i__6].r * b[i__7].r - a[i__6].i *
00480 b[i__7].i, q__2.i = a[i__6].r * b[
00481 i__7].i + a[i__6].i * b[i__7].r;
00482 q__1.r = b[i__5].r - q__2.r, q__1.i = b[i__5]
00483 .i - q__2.i;
00484 b[i__4].r = q__1.r, b[i__4].i = q__1.i;
00485
00486 }
00487 }
00488
00489 }
00490 if (nounit) {
00491 c_div(&q__1, &c_b1, &a[j + j * a_dim1]);
00492 temp.r = q__1.r, temp.i = q__1.i;
00493 i__2 = *m;
00494 for (i__ = 1; i__ <= i__2; ++i__) {
00495 i__3 = i__ + j * b_dim1;
00496 i__4 = i__ + j * b_dim1;
00497 q__1.r = temp.r * b[i__4].r - temp.i * b[i__4].i,
00498 q__1.i = temp.r * b[i__4].i + temp.i * b[
00499 i__4].r;
00500 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00501
00502 }
00503 }
00504
00505 }
00506 } else {
00507 for (j = *n; j >= 1; --j) {
00508 if (alpha->r != 1.f || alpha->i != 0.f) {
00509 i__1 = *m;
00510 for (i__ = 1; i__ <= i__1; ++i__) {
00511 i__2 = i__ + j * b_dim1;
00512 i__3 = i__ + j * b_dim1;
00513 q__1.r = alpha->r * b[i__3].r - alpha->i * b[i__3]
00514 .i, q__1.i = alpha->r * b[i__3].i +
00515 alpha->i * b[i__3].r;
00516 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00517
00518 }
00519 }
00520 i__1 = *n;
00521 for (k = j + 1; k <= i__1; ++k) {
00522 i__2 = k + j * a_dim1;
00523 if (a[i__2].r != 0.f || a[i__2].i != 0.f) {
00524 i__2 = *m;
00525 for (i__ = 1; i__ <= i__2; ++i__) {
00526 i__3 = i__ + j * b_dim1;
00527 i__4 = i__ + j * b_dim1;
00528 i__5 = k + j * a_dim1;
00529 i__6 = i__ + k * b_dim1;
00530 q__2.r = a[i__5].r * b[i__6].r - a[i__5].i *
00531 b[i__6].i, q__2.i = a[i__5].r * b[
00532 i__6].i + a[i__5].i * b[i__6].r;
00533 q__1.r = b[i__4].r - q__2.r, q__1.i = b[i__4]
00534 .i - q__2.i;
00535 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00536
00537 }
00538 }
00539
00540 }
00541 if (nounit) {
00542 c_div(&q__1, &c_b1, &a[j + j * a_dim1]);
00543 temp.r = q__1.r, temp.i = q__1.i;
00544 i__1 = *m;
00545 for (i__ = 1; i__ <= i__1; ++i__) {
00546 i__2 = i__ + j * b_dim1;
00547 i__3 = i__ + j * b_dim1;
00548 q__1.r = temp.r * b[i__3].r - temp.i * b[i__3].i,
00549 q__1.i = temp.r * b[i__3].i + temp.i * b[
00550 i__3].r;
00551 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00552
00553 }
00554 }
00555
00556 }
00557 }
00558 } else {
00559
00560
00561
00562
00563 if (upper) {
00564 for (k = *n; k >= 1; --k) {
00565 if (nounit) {
00566 if (noconj) {
00567 c_div(&q__1, &c_b1, &a[k + k * a_dim1]);
00568 temp.r = q__1.r, temp.i = q__1.i;
00569 } else {
00570 r_cnjg(&q__2, &a[k + k * a_dim1]);
00571 c_div(&q__1, &c_b1, &q__2);
00572 temp.r = q__1.r, temp.i = q__1.i;
00573 }
00574 i__1 = *m;
00575 for (i__ = 1; i__ <= i__1; ++i__) {
00576 i__2 = i__ + k * b_dim1;
00577 i__3 = i__ + k * b_dim1;
00578 q__1.r = temp.r * b[i__3].r - temp.i * b[i__3].i,
00579 q__1.i = temp.r * b[i__3].i + temp.i * b[
00580 i__3].r;
00581 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00582
00583 }
00584 }
00585 i__1 = k - 1;
00586 for (j = 1; j <= i__1; ++j) {
00587 i__2 = j + k * a_dim1;
00588 if (a[i__2].r != 0.f || a[i__2].i != 0.f) {
00589 if (noconj) {
00590 i__2 = j + k * a_dim1;
00591 temp.r = a[i__2].r, temp.i = a[i__2].i;
00592 } else {
00593 r_cnjg(&q__1, &a[j + k * a_dim1]);
00594 temp.r = q__1.r, temp.i = q__1.i;
00595 }
00596 i__2 = *m;
00597 for (i__ = 1; i__ <= i__2; ++i__) {
00598 i__3 = i__ + j * b_dim1;
00599 i__4 = i__ + j * b_dim1;
00600 i__5 = i__ + k * b_dim1;
00601 q__2.r = temp.r * b[i__5].r - temp.i * b[i__5]
00602 .i, q__2.i = temp.r * b[i__5].i +
00603 temp.i * b[i__5].r;
00604 q__1.r = b[i__4].r - q__2.r, q__1.i = b[i__4]
00605 .i - q__2.i;
00606 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00607
00608 }
00609 }
00610
00611 }
00612 if (alpha->r != 1.f || alpha->i != 0.f) {
00613 i__1 = *m;
00614 for (i__ = 1; i__ <= i__1; ++i__) {
00615 i__2 = i__ + k * b_dim1;
00616 i__3 = i__ + k * b_dim1;
00617 q__1.r = alpha->r * b[i__3].r - alpha->i * b[i__3]
00618 .i, q__1.i = alpha->r * b[i__3].i +
00619 alpha->i * b[i__3].r;
00620 b[i__2].r = q__1.r, b[i__2].i = q__1.i;
00621
00622 }
00623 }
00624
00625 }
00626 } else {
00627 i__1 = *n;
00628 for (k = 1; k <= i__1; ++k) {
00629 if (nounit) {
00630 if (noconj) {
00631 c_div(&q__1, &c_b1, &a[k + k * a_dim1]);
00632 temp.r = q__1.r, temp.i = q__1.i;
00633 } else {
00634 r_cnjg(&q__2, &a[k + k * a_dim1]);
00635 c_div(&q__1, &c_b1, &q__2);
00636 temp.r = q__1.r, temp.i = q__1.i;
00637 }
00638 i__2 = *m;
00639 for (i__ = 1; i__ <= i__2; ++i__) {
00640 i__3 = i__ + k * b_dim1;
00641 i__4 = i__ + k * b_dim1;
00642 q__1.r = temp.r * b[i__4].r - temp.i * b[i__4].i,
00643 q__1.i = temp.r * b[i__4].i + temp.i * b[
00644 i__4].r;
00645 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00646
00647 }
00648 }
00649 i__2 = *n;
00650 for (j = k + 1; j <= i__2; ++j) {
00651 i__3 = j + k * a_dim1;
00652 if (a[i__3].r != 0.f || a[i__3].i != 0.f) {
00653 if (noconj) {
00654 i__3 = j + k * a_dim1;
00655 temp.r = a[i__3].r, temp.i = a[i__3].i;
00656 } else {
00657 r_cnjg(&q__1, &a[j + k * a_dim1]);
00658 temp.r = q__1.r, temp.i = q__1.i;
00659 }
00660 i__3 = *m;
00661 for (i__ = 1; i__ <= i__3; ++i__) {
00662 i__4 = i__ + j * b_dim1;
00663 i__5 = i__ + j * b_dim1;
00664 i__6 = i__ + k * b_dim1;
00665 q__2.r = temp.r * b[i__6].r - temp.i * b[i__6]
00666 .i, q__2.i = temp.r * b[i__6].i +
00667 temp.i * b[i__6].r;
00668 q__1.r = b[i__5].r - q__2.r, q__1.i = b[i__5]
00669 .i - q__2.i;
00670 b[i__4].r = q__1.r, b[i__4].i = q__1.i;
00671
00672 }
00673 }
00674
00675 }
00676 if (alpha->r != 1.f || alpha->i != 0.f) {
00677 i__2 = *m;
00678 for (i__ = 1; i__ <= i__2; ++i__) {
00679 i__3 = i__ + k * b_dim1;
00680 i__4 = i__ + k * b_dim1;
00681 q__1.r = alpha->r * b[i__4].r - alpha->i * b[i__4]
00682 .i, q__1.i = alpha->r * b[i__4].i +
00683 alpha->i * b[i__4].r;
00684 b[i__3].r = q__1.r, b[i__3].i = q__1.i;
00685
00686 }
00687 }
00688
00689 }
00690 }
00691 }
00692 }
00693
00694 return 0;
00695
00696
00697
00698 }