source: orange/source/orange/linpack/dtrsl.c @ 8069:5e19a5dcaac0

Revision 8069:5e19a5dcaac0, 4.9 KB checked in by ales_erjavec <ales.erjavec@…>, 3 years ago (diff)
Line 
1/* dtrsl.f -- translated by f2c (version 20090411).
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
10        http://www.netlib.org/f2c/libf2c.zip
11*/
12
13//#include "f2c"
14#include "linpack.h"
15
16/* Table of constant values */
17
18static integer c__1 = 1;
19
20/* Subroutine */ int dtrsl_(doublereal *t, integer *ldt, integer *n, 
21    doublereal *b, integer *job, integer *info)
22{
23    /* System generated locals */
24    integer t_dim1, t_offset, i__1, i__2;
25
26    /* Local variables */
27    static integer j, jj, case__;
28    extern doublereal ddot_(integer *, doublereal *, integer *, doublereal *, 
29        integer *);
30    static doublereal temp;
31    extern /* Subroutine */ int daxpy_(integer *, doublereal *, doublereal *, 
32        integer *, doublereal *, integer *);
33
34
35
36/*     dtrsl solves systems of the form */
37
38/*                   t * x = b */
39/*     or */
40/*                   trans(t) * x = b */
41
42/*     where t is a triangular matrix of order n. here trans(t) */
43/*     denotes the transpose of the matrix t. */
44
45/*     on entry */
46
47/*         t         double precision(ldt,n) */
48/*                   t contains the matrix of the system. the zero */
49/*                   elements of the matrix are not referenced, and */
50/*                   the corresponding elements of the array can be */
51/*                   used to store other information. */
52
53/*         ldt       integer */
54/*                   ldt is the leading dimension of the array t. */
55
56/*         n         integer */
57/*                   n is the order of the system. */
58
59/*         b         double precision(n). */
60/*                   b contains the right hand side of the system. */
61
62/*         job       integer */
63/*                   job specifies what kind of system is to be solved. */
64/*                   if job is */
65
66/*                        00   solve t*x=b, t lower triangular, */
67/*                        01   solve t*x=b, t upper triangular, */
68/*                        10   solve trans(t)*x=b, t lower triangular, */
69/*                        11   solve trans(t)*x=b, t upper triangular. */
70
71/*     on return */
72
73/*         b         b contains the solution, if info .eq. 0. */
74/*                   otherwise b is unaltered. */
75
76/*         info      integer */
77/*                   info contains zero if the system is nonsingular. */
78/*                   otherwise info contains the index of */
79/*                   the first zero diagonal element of t. */
80
81/*     linpack. this version dated 08/14/78 . */
82/*     g. w. stewart, university of maryland, argonne national lab. */
83
84/*     subroutines and functions */
85
86/*     blas daxpy,ddot */
87/*     fortran mod */
88
89/*     internal variables */
90
91
92/*     begin block permitting ...exits to 150 */
93
94/*        check for zero diagonal elements. */
95
96    /* Parameter adjustments */
97    t_dim1 = *ldt;
98    t_offset = 1 + t_dim1;
99    t -= t_offset;
100    --b;
101
102    /* Function Body */
103    i__1 = *n;
104    for (*info = 1; *info <= i__1; ++(*info)) {
105/*     ......exit */
106    if (t[*info + *info * t_dim1] == 0.) {
107        goto L150;
108    }
109/* L10: */
110    }
111    *info = 0;
112
113/*        determine the task and go to it. */
114
115    case__ = 1;
116    if (*job % 10 != 0) {
117    case__ = 2;
118    }
119    if (*job % 100 / 10 != 0) {
120    case__ += 2;
121    }
122    switch (case__) {
123    case 1goto L20;
124    case 2goto L50;
125    case 3goto L80;
126    case 4goto L110;
127    }
128
129/*        solve t*x=b for t lower triangular */
130
131L20:
132    b[1] /= t[t_dim1 + 1];
133    if (*n < 2) {
134    goto L40;
135    }
136    i__1 = *n;
137    for (j = 2; j <= i__1; ++j) {
138    temp = -b[j - 1];
139    i__2 = *n - j + 1;
140    daxpy_(&i__2, &temp, &t[j + (j - 1) * t_dim1], &c__1, &b[j], &c__1);
141    b[j] /= t[j + j * t_dim1];
142/* L30: */
143    }
144L40:
145    goto L140;
146
147/*        solve t*x=b for t upper triangular. */
148
149L50:
150    b[*n] /= t[*n + *n * t_dim1];
151    if (*n < 2) {
152    goto L70;
153    }
154    i__1 = *n;
155    for (jj = 2; jj <= i__1; ++jj) {
156    j = *n - jj + 1;
157    temp = -b[j + 1];
158    daxpy_(&j, &temp, &t[(j + 1) * t_dim1 + 1], &c__1, &b[1], &c__1);
159    b[j] /= t[j + j * t_dim1];
160/* L60: */
161    }
162L70:
163    goto L140;
164
165/*        solve trans(t)*x=b for t lower triangular. */
166
167L80:
168    b[*n] /= t[*n + *n * t_dim1];
169    if (*n < 2) {
170    goto L100;
171    }
172    i__1 = *n;
173    for (jj = 2; jj <= i__1; ++jj) {
174    j = *n - jj + 1;
175    i__2 = jj - 1;
176    b[j] -= ddot_(&i__2, &t[j + 1 + j * t_dim1], &c__1, &b[j + 1], &c__1);
177    b[j] /= t[j + j * t_dim1];
178/* L90: */
179    }
180L100:
181    goto L140;
182
183/*        solve trans(t)*x=b for t upper triangular. */
184
185L110:
186    b[1] /= t[t_dim1 + 1];
187    if (*n < 2) {
188    goto L130;
189    }
190    i__1 = *n;
191    for (j = 2; j <= i__1; ++j) {
192    i__2 = j - 1;
193    b[j] -= ddot_(&i__2, &t[j * t_dim1 + 1], &c__1, &b[1], &c__1);
194    b[j] /= t[j + j * t_dim1];
195/* L120: */
196    }
197L130:
198L140:
199L150:
200    return 0;
201} /* dtrsl_ */
202
Note: See TracBrowser for help on using the repository browser.