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