37 #ifndef TEMPLATE_LAPACK_LATRS_HEADER 38 #define TEMPLATE_LAPACK_LATRS_HEADER 43 normin,
const integer *n,
const Treal *a,
const integer *lda, Treal *x,
44 Treal *scale, Treal *cnorm,
integer *info)
213 integer a_dim1, a_offset, i__1, i__2, i__3;
214 Treal d__1, d__2, d__3;
219 Treal tmax, tjjs, xmax, grow, sumj;
231 #define a_ref(a_1,a_2) a[(a_2)*a_dim1 + a_1] 235 a_offset = 1 + a_dim1 * 1;
278 bignum = 1. / smlnum;
290 for (j = 1; j <= i__1; ++j) {
300 for (j = 1; j <= i__1; ++j) {
314 if (tmax <= bignum) {
317 tscal = 1. / (smlnum * tmax);
318 dscal_(n, &tscal, &cnorm[1], &c__1);
325 xmax = (d__1 = x[j],
absMACRO(d__1));
357 for (j = jfirst; i__2 < 0 ? j >= i__1 : j <= i__1; j += i__2) {
361 if (grow <= smlnum) {
369 d__1 = xbnd, d__2 =
minMACRO(1.,tjj) * grow;
371 if (tjj + cnorm[j] >= smlnum) {
375 grow *= tjj / (tjj + cnorm[j]);
392 d__1 = 1., d__2 = 1. /
maxMACRO(xbnd,smlnum);
396 for (j = jfirst; i__1 < 0 ? j >= i__2 : j <= i__2; j += i__1) {
400 if (grow <= smlnum) {
406 grow *= 1. / (cnorm[j] + 1.);
443 for (j = jfirst; i__2 < 0 ? j >= i__1 : j <= i__1; j += i__2) {
447 if (grow <= smlnum) {
455 d__1 = grow, d__2 = xbnd / xj;
474 d__1 = 1., d__2 = 1. /
maxMACRO(xbnd,smlnum);
478 for (j = jfirst; i__1 < 0 ? j >= i__2 : j <= i__2; j += i__1) {
482 if (grow <= smlnum) {
497 if (grow * tscal > smlnum) {
512 *scale = bignum / xmax;
513 dscal_(n, scale, &x[1], &c__1);
523 for (j = jfirst; i__2 < 0 ? j >= i__1 : j <= i__1; j += i__2) {
529 tjjs =
a_ref(j, j) * tscal;
542 if (xj > tjj * bignum) {
547 dscal_(n, &rec, &x[1], &c__1);
554 }
else if (tjj > 0.) {
558 if (xj > tjj * bignum) {
563 rec = tjj * bignum / xj;
571 dscal_(n, &rec, &x[1], &c__1);
583 for (i__ = 1; i__ <= i__3; ++i__) {
599 if (cnorm[j] > (bignum - xmax) * rec) {
604 dscal_(n, &rec, &x[1], &c__1);
607 }
else if (xj * cnorm[j] > bignum - xmax) {
611 dscal_(n, &c_b36, &x[1], &c__1);
622 d__1 = -x[j] * tscal;
627 xmax = (d__1 = x[i__],
absMACRO(d__1));
636 d__1 = -x[j] * tscal;
637 daxpy_(&i__3, &d__1, &
a_ref(j + 1, j), &c__1, &x[j +
641 xmax = (d__1 = x[i__],
absMACRO(d__1));
653 for (j = jfirst; i__1 < 0 ? j >= i__2 : j <= i__2; j += i__1) {
661 if (cnorm[j] > (bignum - xj) * rec) {
667 tjjs =
a_ref(j, j) * tscal;
677 d__1 = 1., d__2 = rec * tjj;
682 dscal_(n, &rec, &x[1], &c__1);
696 sumj =
ddot_(&i__3, &
a_ref(1, j), &c__1, &x[1], &c__1)
700 sumj =
ddot_(&i__3, &
a_ref(j + 1, j), &c__1, &x[j + 1]
709 for (i__ = 1; i__ <= i__3; ++i__) {
710 sumj +=
a_ref(i__, j) * uscal * x[i__];
715 for (i__ = j + 1; i__ <= i__3; ++i__) {
716 sumj +=
a_ref(i__, j) * uscal * x[i__];
722 if (uscal == tscal) {
730 tjjs =
a_ref(j, j) * tscal;
746 if (xj > tjj * bignum) {
751 dscal_(n, &rec, &x[1], &c__1);
757 }
else if (tjj > 0.) {
761 if (xj > tjj * bignum) {
765 rec = tjj * bignum / xj;
766 dscal_(n, &rec, &x[1], &c__1);
777 for (i__ = 1; i__ <= i__3; ++i__) {
792 x[j] = x[j] / tjjs - sumj;
795 d__2 = xmax, d__3 = (d__1 = x[j],
absMACRO(d__1));
807 dscal_(n, &d__1, &cnorm[1], &c__1);
#define absMACRO(x)
Definition: template_blas_common.h:47
integer template_blas_idamax(const integer *n, const Treal *dx, const integer *incx)
Definition: template_blas_idamax.h:42
int integer
Definition: template_blas_common.h:40
#define maxMACRO(a, b)
Definition: template_blas_common.h:45
void daxpy_(const int *n, const double *da, const double *dx, const int *incx, double *dy, const int *incy)
Treal template_blas_asum(const integer *n, const Treal *dx, const integer *incx)
Definition: template_blas_asum.h:42
#define minMACRO(a, b)
Definition: template_blas_common.h:46
int template_blas_erbla(const char *srname, integer *info)
Definition: template_blas_common.cc:146
int template_lapack_latrs(const char *uplo, const char *trans, const char *diag, const char *normin, const integer *n, const Treal *a, const integer *lda, Treal *x, Treal *scale, Treal *cnorm, integer *info)
Definition: template_lapack_latrs.h:42
void dscal_(const int *n, const double *da, double *dx, const int *incx)
int template_blas_trsv(const char *uplo, const char *trans, const char *diag, const integer *n, const Treal *a, const integer *lda, Treal *x, const integer *incx)
Definition: template_blas_trsv.h:42
Treal template_lapack_lamch(const char *cmach, Treal dummyReal)
Definition: template_lapack_lamch.h:202
bool logical
Definition: template_blas_common.h:41
double ddot_(const int *n, const double *dx, const int *incx, const double *dy, const int *incy)
logical template_blas_lsame(const char *ca, const char *cb)
Definition: template_blas_common.cc:46