00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035 #ifndef TEMPLATE_LAPACK_LARFB_HEADER
00036 #define TEMPLATE_LAPACK_LARFB_HEADER
00037
00038
00039 template<class Treal>
00040 int template_lapack_larfb(const char *side, const char *trans, const char *direct, const char *
00041 storev, const integer *m, const integer *n, const integer *k, const Treal *v, const integer *
00042 ldv, const Treal *t, const integer *ldt, Treal *c__, const integer *ldc,
00043 Treal *work, const integer *ldwork)
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 integer c__1 = 1;
00131 Treal c_b14 = 1.;
00132 Treal c_b25 = -1.;
00133
00134
00135 integer c_dim1, c_offset, t_dim1, t_offset, v_dim1, v_offset, work_dim1,
00136 work_offset, i__1, i__2;
00137
00138 integer i__, j;
00139 char transt[1];
00140 #define work_ref(a_1,a_2) work[(a_2)*work_dim1 + a_1]
00141 #define c___ref(a_1,a_2) c__[(a_2)*c_dim1 + a_1]
00142 #define v_ref(a_1,a_2) v[(a_2)*v_dim1 + a_1]
00143
00144
00145 v_dim1 = *ldv;
00146 v_offset = 1 + v_dim1 * 1;
00147 v -= v_offset;
00148 t_dim1 = *ldt;
00149 t_offset = 1 + t_dim1 * 1;
00150 t -= t_offset;
00151 c_dim1 = *ldc;
00152 c_offset = 1 + c_dim1 * 1;
00153 c__ -= c_offset;
00154 work_dim1 = *ldwork;
00155 work_offset = 1 + work_dim1 * 1;
00156 work -= work_offset;
00157
00158
00159 if (*m <= 0 || *n <= 0) {
00160 return 0;
00161 }
00162
00163 if (template_blas_lsame(trans, "N")) {
00164 *(unsigned char *)transt = 'T';
00165 } else {
00166 *(unsigned char *)transt = 'N';
00167 }
00168
00169 if (template_blas_lsame(storev, "C")) {
00170
00171 if (template_blas_lsame(direct, "F")) {
00172
00173
00174
00175
00176
00177 if (template_blas_lsame(side, "L")) {
00178
00179
00180
00181
00182
00183
00184
00185
00186 i__1 = *k;
00187 for (j = 1; j <= i__1; ++j) {
00188 template_blas_copy(n, &c___ref(j, 1), ldc, &work_ref(1, j), &c__1);
00189
00190 }
00191
00192
00193
00194 template_blas_trmm("Right", "Lower", "No transpose", "Unit", n, k, &c_b14,
00195 &v[v_offset], ldv, &work[work_offset], ldwork);
00196 if (*m > *k) {
00197
00198
00199
00200 i__1 = *m - *k;
00201 template_blas_gemm("Transpose", "No transpose", n, k, &i__1, &c_b14, &
00202 c___ref(*k + 1, 1), ldc, &v_ref(*k + 1, 1), ldv, &
00203 c_b14, &work[work_offset], ldwork);
00204 }
00205
00206
00207
00208 template_blas_trmm("Right", "Upper", transt, "Non-unit", n, k, &c_b14, &t[
00209 t_offset], ldt, &work[work_offset], ldwork);
00210
00211
00212
00213 if (*m > *k) {
00214
00215
00216
00217 i__1 = *m - *k;
00218 template_blas_gemm("No transpose", "Transpose", &i__1, n, k, &c_b25, &
00219 v_ref(*k + 1, 1), ldv, &work[work_offset], ldwork,
00220 &c_b14, &c___ref(*k + 1, 1), ldc);
00221 }
00222
00223
00224
00225 template_blas_trmm("Right", "Lower", "Transpose", "Unit", n, k, &c_b14, &
00226 v[v_offset], ldv, &work[work_offset], ldwork);
00227
00228
00229
00230 i__1 = *k;
00231 for (j = 1; j <= i__1; ++j) {
00232 i__2 = *n;
00233 for (i__ = 1; i__ <= i__2; ++i__) {
00234 c___ref(j, i__) = c___ref(j, i__) - work_ref(i__, j);
00235
00236 }
00237
00238 }
00239
00240 } else if (template_blas_lsame(side, "R")) {
00241
00242
00243
00244
00245
00246
00247
00248 i__1 = *k;
00249 for (j = 1; j <= i__1; ++j) {
00250 template_blas_copy(m, &c___ref(1, j), &c__1, &work_ref(1, j), &c__1);
00251
00252 }
00253
00254
00255
00256 template_blas_trmm("Right", "Lower", "No transpose", "Unit", m, k, &c_b14,
00257 &v[v_offset], ldv, &work[work_offset], ldwork);
00258 if (*n > *k) {
00259
00260
00261
00262 i__1 = *n - *k;
00263 template_blas_gemm("No transpose", "No transpose", m, k, &i__1, &
00264 c_b14, &c___ref(1, *k + 1), ldc, &v_ref(*k + 1, 1)
00265 , ldv, &c_b14, &work[work_offset], ldwork);
00266 }
00267
00268
00269
00270 template_blas_trmm("Right", "Upper", trans, "Non-unit", m, k, &c_b14, &t[
00271 t_offset], ldt, &work[work_offset], ldwork);
00272
00273
00274
00275 if (*n > *k) {
00276
00277
00278
00279 i__1 = *n - *k;
00280 template_blas_gemm("No transpose", "Transpose", m, &i__1, k, &c_b25, &
00281 work[work_offset], ldwork, &v_ref(*k + 1, 1), ldv,
00282 &c_b14, &c___ref(1, *k + 1), ldc);
00283 }
00284
00285
00286
00287 template_blas_trmm("Right", "Lower", "Transpose", "Unit", m, k, &c_b14, &
00288 v[v_offset], ldv, &work[work_offset], ldwork);
00289
00290
00291
00292 i__1 = *k;
00293 for (j = 1; j <= i__1; ++j) {
00294 i__2 = *m;
00295 for (i__ = 1; i__ <= i__2; ++i__) {
00296 c___ref(i__, j) = c___ref(i__, j) - work_ref(i__, j);
00297
00298 }
00299
00300 }
00301 }
00302
00303 } else {
00304
00305
00306
00307
00308
00309 if (template_blas_lsame(side, "L")) {
00310
00311
00312
00313
00314
00315
00316
00317
00318 i__1 = *k;
00319 for (j = 1; j <= i__1; ++j) {
00320 template_blas_copy(n, &c___ref(*m - *k + j, 1), ldc, &work_ref(1, j),
00321 &c__1);
00322
00323 }
00324
00325
00326
00327 template_blas_trmm("Right", "Upper", "No transpose", "Unit", n, k, &c_b14,
00328 &v_ref(*m - *k + 1, 1), ldv, &work[work_offset],
00329 ldwork);
00330 if (*m > *k) {
00331
00332
00333
00334 i__1 = *m - *k;
00335 template_blas_gemm("Transpose", "No transpose", n, k, &i__1, &c_b14, &
00336 c__[c_offset], ldc, &v[v_offset], ldv, &c_b14, &
00337 work[work_offset], ldwork);
00338 }
00339
00340
00341
00342 template_blas_trmm("Right", "Lower", transt, "Non-unit", n, k, &c_b14, &t[
00343 t_offset], ldt, &work[work_offset], ldwork);
00344
00345
00346
00347 if (*m > *k) {
00348
00349
00350
00351 i__1 = *m - *k;
00352 template_blas_gemm("No transpose", "Transpose", &i__1, n, k, &c_b25, &
00353 v[v_offset], ldv, &work[work_offset], ldwork, &
00354 c_b14, &c__[c_offset], ldc)
00355 ;
00356 }
00357
00358
00359
00360 template_blas_trmm("Right", "Upper", "Transpose", "Unit", n, k, &c_b14, &
00361 v_ref(*m - *k + 1, 1), ldv, &work[work_offset],
00362 ldwork);
00363
00364
00365
00366 i__1 = *k;
00367 for (j = 1; j <= i__1; ++j) {
00368 i__2 = *n;
00369 for (i__ = 1; i__ <= i__2; ++i__) {
00370 c___ref(*m - *k + j, i__) = c___ref(*m - *k + j, i__)
00371 - work_ref(i__, j);
00372
00373 }
00374
00375 }
00376
00377 } else if (template_blas_lsame(side, "R")) {
00378
00379
00380
00381
00382
00383
00384
00385 i__1 = *k;
00386 for (j = 1; j <= i__1; ++j) {
00387 template_blas_copy(m, &c___ref(1, *n - *k + j), &c__1, &work_ref(1, j)
00388 , &c__1);
00389
00390 }
00391
00392
00393
00394 template_blas_trmm("Right", "Upper", "No transpose", "Unit", m, k, &c_b14,
00395 &v_ref(*n - *k + 1, 1), ldv, &work[work_offset],
00396 ldwork);
00397 if (*n > *k) {
00398
00399
00400
00401 i__1 = *n - *k;
00402 template_blas_gemm("No transpose", "No transpose", m, k, &i__1, &
00403 c_b14, &c__[c_offset], ldc, &v[v_offset], ldv, &
00404 c_b14, &work[work_offset], ldwork);
00405 }
00406
00407
00408
00409 template_blas_trmm("Right", "Lower", trans, "Non-unit", m, k, &c_b14, &t[
00410 t_offset], ldt, &work[work_offset], ldwork);
00411
00412
00413
00414 if (*n > *k) {
00415
00416
00417
00418 i__1 = *n - *k;
00419 template_blas_gemm("No transpose", "Transpose", m, &i__1, k, &c_b25, &
00420 work[work_offset], ldwork, &v[v_offset], ldv, &
00421 c_b14, &c__[c_offset], ldc)
00422 ;
00423 }
00424
00425
00426
00427 template_blas_trmm("Right", "Upper", "Transpose", "Unit", m, k, &c_b14, &
00428 v_ref(*n - *k + 1, 1), ldv, &work[work_offset],
00429 ldwork);
00430
00431
00432
00433 i__1 = *k;
00434 for (j = 1; j <= i__1; ++j) {
00435 i__2 = *m;
00436 for (i__ = 1; i__ <= i__2; ++i__) {
00437 c___ref(i__, *n - *k + j) = c___ref(i__, *n - *k + j)
00438 - work_ref(i__, j);
00439
00440 }
00441
00442 }
00443 }
00444 }
00445
00446 } else if (template_blas_lsame(storev, "R")) {
00447
00448 if (template_blas_lsame(direct, "F")) {
00449
00450
00451
00452
00453 if (template_blas_lsame(side, "L")) {
00454
00455
00456
00457
00458
00459
00460
00461
00462 i__1 = *k;
00463 for (j = 1; j <= i__1; ++j) {
00464 template_blas_copy(n, &c___ref(j, 1), ldc, &work_ref(1, j), &c__1);
00465
00466 }
00467
00468
00469
00470 template_blas_trmm("Right", "Upper", "Transpose", "Unit", n, k, &c_b14, &
00471 v[v_offset], ldv, &work[work_offset], ldwork);
00472 if (*m > *k) {
00473
00474
00475
00476 i__1 = *m - *k;
00477 template_blas_gemm("Transpose", "Transpose", n, k, &i__1, &c_b14, &
00478 c___ref(*k + 1, 1), ldc, &v_ref(1, *k + 1), ldv, &
00479 c_b14, &work[work_offset], ldwork);
00480 }
00481
00482
00483
00484 template_blas_trmm("Right", "Upper", transt, "Non-unit", n, k, &c_b14, &t[
00485 t_offset], ldt, &work[work_offset], ldwork);
00486
00487
00488
00489 if (*m > *k) {
00490
00491
00492
00493 i__1 = *m - *k;
00494 template_blas_gemm("Transpose", "Transpose", &i__1, n, k, &c_b25, &
00495 v_ref(1, *k + 1), ldv, &work[work_offset], ldwork,
00496 &c_b14, &c___ref(*k + 1, 1), ldc);
00497 }
00498
00499
00500
00501 template_blas_trmm("Right", "Upper", "No transpose", "Unit", n, k, &c_b14,
00502 &v[v_offset], ldv, &work[work_offset], ldwork);
00503
00504
00505
00506 i__1 = *k;
00507 for (j = 1; j <= i__1; ++j) {
00508 i__2 = *n;
00509 for (i__ = 1; i__ <= i__2; ++i__) {
00510 c___ref(j, i__) = c___ref(j, i__) - work_ref(i__, j);
00511
00512 }
00513
00514 }
00515
00516 } else if (template_blas_lsame(side, "R")) {
00517
00518
00519
00520
00521
00522
00523
00524 i__1 = *k;
00525 for (j = 1; j <= i__1; ++j) {
00526 template_blas_copy(m, &c___ref(1, j), &c__1, &work_ref(1, j), &c__1);
00527
00528 }
00529
00530
00531
00532 template_blas_trmm("Right", "Upper", "Transpose", "Unit", m, k, &c_b14, &
00533 v[v_offset], ldv, &work[work_offset], ldwork);
00534 if (*n > *k) {
00535
00536
00537
00538 i__1 = *n - *k;
00539 template_blas_gemm("No transpose", "Transpose", m, k, &i__1, &c_b14, &
00540 c___ref(1, *k + 1), ldc, &v_ref(1, *k + 1), ldv, &
00541 c_b14, &work[work_offset], ldwork);
00542 }
00543
00544
00545
00546 template_blas_trmm("Right", "Upper", trans, "Non-unit", m, k, &c_b14, &t[
00547 t_offset], ldt, &work[work_offset], ldwork);
00548
00549
00550
00551 if (*n > *k) {
00552
00553
00554
00555 i__1 = *n - *k;
00556 template_blas_gemm("No transpose", "No transpose", m, &i__1, k, &
00557 c_b25, &work[work_offset], ldwork, &v_ref(1, *k +
00558 1), ldv, &c_b14, &c___ref(1, *k + 1), ldc);
00559 }
00560
00561
00562
00563 template_blas_trmm("Right", "Upper", "No transpose", "Unit", m, k, &c_b14,
00564 &v[v_offset], ldv, &work[work_offset], ldwork);
00565
00566
00567
00568 i__1 = *k;
00569 for (j = 1; j <= i__1; ++j) {
00570 i__2 = *m;
00571 for (i__ = 1; i__ <= i__2; ++i__) {
00572 c___ref(i__, j) = c___ref(i__, j) - work_ref(i__, j);
00573
00574 }
00575
00576 }
00577
00578 }
00579
00580 } else {
00581
00582
00583
00584
00585 if (template_blas_lsame(side, "L")) {
00586
00587
00588
00589
00590
00591
00592
00593
00594 i__1 = *k;
00595 for (j = 1; j <= i__1; ++j) {
00596 template_blas_copy(n, &c___ref(*m - *k + j, 1), ldc, &work_ref(1, j),
00597 &c__1);
00598
00599 }
00600
00601
00602
00603 template_blas_trmm("Right", "Lower", "Transpose", "Unit", n, k, &c_b14, &
00604 v_ref(1, *m - *k + 1), ldv, &work[work_offset],
00605 ldwork);
00606 if (*m > *k) {
00607
00608
00609
00610 i__1 = *m - *k;
00611 template_blas_gemm("Transpose", "Transpose", n, k, &i__1, &c_b14, &
00612 c__[c_offset], ldc, &v[v_offset], ldv, &c_b14, &
00613 work[work_offset], ldwork);
00614 }
00615
00616
00617
00618 template_blas_trmm("Right", "Lower", transt, "Non-unit", n, k, &c_b14, &t[
00619 t_offset], ldt, &work[work_offset], ldwork);
00620
00621
00622
00623 if (*m > *k) {
00624
00625
00626
00627 i__1 = *m - *k;
00628 template_blas_gemm("Transpose", "Transpose", &i__1, n, k, &c_b25, &v[
00629 v_offset], ldv, &work[work_offset], ldwork, &
00630 c_b14, &c__[c_offset], ldc);
00631 }
00632
00633
00634
00635 template_blas_trmm("Right", "Lower", "No transpose", "Unit", n, k, &c_b14,
00636 &v_ref(1, *m - *k + 1), ldv, &work[work_offset],
00637 ldwork);
00638
00639
00640
00641 i__1 = *k;
00642 for (j = 1; j <= i__1; ++j) {
00643 i__2 = *n;
00644 for (i__ = 1; i__ <= i__2; ++i__) {
00645 c___ref(*m - *k + j, i__) = c___ref(*m - *k + j, i__)
00646 - work_ref(i__, j);
00647
00648 }
00649
00650 }
00651
00652 } else if (template_blas_lsame(side, "R")) {
00653
00654
00655
00656
00657
00658
00659
00660 i__1 = *k;
00661 for (j = 1; j <= i__1; ++j) {
00662 template_blas_copy(m, &c___ref(1, *n - *k + j), &c__1, &work_ref(1, j)
00663 , &c__1);
00664
00665 }
00666
00667
00668
00669 template_blas_trmm("Right", "Lower", "Transpose", "Unit", m, k, &c_b14, &
00670 v_ref(1, *n - *k + 1), ldv, &work[work_offset],
00671 ldwork);
00672 if (*n > *k) {
00673
00674
00675
00676 i__1 = *n - *k;
00677 template_blas_gemm("No transpose", "Transpose", m, k, &i__1, &c_b14, &
00678 c__[c_offset], ldc, &v[v_offset], ldv, &c_b14, &
00679 work[work_offset], ldwork);
00680 }
00681
00682
00683
00684 template_blas_trmm("Right", "Lower", trans, "Non-unit", m, k, &c_b14, &t[
00685 t_offset], ldt, &work[work_offset], ldwork);
00686
00687
00688
00689 if (*n > *k) {
00690
00691
00692
00693 i__1 = *n - *k;
00694 template_blas_gemm("No transpose", "No transpose", m, &i__1, k, &
00695 c_b25, &work[work_offset], ldwork, &v[v_offset],
00696 ldv, &c_b14, &c__[c_offset], ldc);
00697 }
00698
00699
00700
00701 template_blas_trmm("Right", "Lower", "No transpose", "Unit", m, k, &c_b14,
00702 &v_ref(1, *n - *k + 1), ldv, &work[work_offset],
00703 ldwork);
00704
00705
00706
00707 i__1 = *k;
00708 for (j = 1; j <= i__1; ++j) {
00709 i__2 = *m;
00710 for (i__ = 1; i__ <= i__2; ++i__) {
00711 c___ref(i__, *n - *k + j) = c___ref(i__, *n - *k + j)
00712 - work_ref(i__, j);
00713
00714 }
00715
00716 }
00717
00718 }
00719
00720 }
00721 }
00722
00723 return 0;
00724
00725
00726
00727 }
00728
00729 #undef v_ref
00730 #undef c___ref
00731 #undef work_ref
00732
00733
00734 #endif