Партнёры:
Хостинг:
NAME
zqrsl - solve the linear system Ax = b for a general matrix
A, which has been QR- factored by xQRDC, and vectors b and
x.
SYNOPSIS
SUBROUTINE DQRSL (DA, LDA, N, K, DQRAUX, DY, DQY, DQTY, DB,
DRESID, DAB, JOB, INFO)
SUBROUTINE SQRSL (SA, LDA, N, K, SQRAUX, SY, SQY, SQTY, SB,
SRESID, SAB, JOB, INFO)
SUBROUTINE ZQRSL (ZA, LDA, N, K, ZQRAUX, CY, CQY, ZQTY, ZB,
ZRESID, ZAB, JOB, INFO)
SUBROUTINE CQRSL (CA, LDA, N, K, CQRAUX, ZY, ZQY, CQTY, CB,
CRESID, CAB, JOB, INFO)
#include <sunperf.h>
void dqrsl(double *dx, int ldx, int n, int k, double *qraux,
double *dy, double *qy, double *qty, double *db,
double *rsd, double *xb, int job, int *info) ;
void sqrsl(float *sx, int ldx, int n, int k, float *qraux,
float *sy, float *qy, float *qty, float *sb, float
*rsd, float *xb, int job, int *info) ;
void zqrsl(doublecomplex *zx, int ldx, int n, int k, doub-
lecomplex *qraux, doublecomplex *zy, doublecomplex
*qy, doublecomplex *qty, doublecomplex *zb, doub-
lecomplex *rsd, doublecomplex *xb, int job, int
*info) ;
void cqrsl(complex *cx, int ldx, int n, int k, complex
*qraux, complex *cy, complex *qy, complex *qty,
complex *b, complex *rsd, complex *xb, int job,
int *info) ;
ARGUMENTS
xA Part of the QR factorization of matrix A as com-
puted by xSRDC.
LDA Leading dimension of the array A as specified in a
dimension or type statement. LDA >= max(1,N).
N Number of rows in the matrix AK where AK is
described below. N >= 0.
K Number of columns in the matrix AK where AK is
described below. K >= 0.
xQRAUX Auxiliary output from xQRDC.
xY Vector to be manipulated by xQRSL.
xQY On exit, QY contains Q * Y if its computation has
been requested in JOB; QY is not referenced if its
computation is not requested.
xQTY On exit, QTY contains QT * Y if its computation
has been requested in JOB; QTY is not referenced
if its computation is not requested.
xB On entry, the right-hand side vector b. On exit,
the solution vector x. B is not referenced if its
computation is not requested.
xRESID On exit, RESID contains the least squares residual
y - AK * b if its computation has been requested.
RESID also is the orthogonal projection of y onto
the orthogonal complement of the column space of
AK. RESID is not referenced if its computation is
not requested.
xAB On exit, AB contains the least squares approxima-
tion AK * b if its computation has been requested.
AB is the orthogonal projection of y onto the
column space of x. AB is not referenced if its
computation is not requested.
JOB Integer in the form abcde; determines which opera-
tion or operations the subroutine will perform:
a <> 0 compute QY
b, c, d, or e <> 0 compute QTY
c <> 0 compute B
d <> 0 compute RESID
e <> 0 compute AB
INFO On exit:
INFO = 0 Subroutine completed normally.
INFO * 0 Returns a value k if the computation of
B has been requested and R is singular; the value
of k is then the index of the first zero element
of R. The matrix AK is constructed from the fac-
tored orthogonal matrix Q and upper triangular
matrix R from xQRDC.
SAMPLE PROGRAM
PROGRAM TEST
IMPLICIT NONE
C
INTEGER IDOB, IDORSD, IDOXB, LDA, N, NCOLA, NOPIV, NROWA
PARAMETER (IDOB = 100)
PARAMETER (IDORSD = 10)
PARAMETER (IDOXB = 1)
PARAMETER (N = 3)
PARAMETER (LDA = N)
PARAMETER (NCOLA = 2)
PARAMETER (NOPIV = 0)
PARAMETER (NROWA = N)
C
DOUBLE PRECISION A(LDA,NCOLA), B(NCOLA), NULL(N), QRAUX(N)
DOUBLE PRECISION RESID(N), WORK(N), Y(N)
INTEGER ICOL, INFO, IROW, JOB, JPIVOT
C
EXTERNAL DQRDC, DQRSL
C
C Initialize the array A to store the matrix A shown below.
C Initialize the array Y to store the vector y shown below.
C
C 1 1 1
C A = 1 0 y = 0
C 0 1 -5
C
DATA A / 1.0D0, 1.0D0, 0.0D0, 1.0D0, 0.0D0, 1.0D0 /
DATA Y / 1.0D0, 0.0D0, -5.0D0 /
C
PRINT 1000
PRINT 1010, ((A(IROW,ICOL), ICOL = 1, NCOLA), IROW = 1, NROWA)
PRINT 1020
PRINT 1030, Y
JOB = NOPIV
CALL DQRDC (A, LDA, NROWA, NCOLA, QRAUX, JPIVOT, WORK, JOB)
JOB = IDOB + IDORSD + IDOXB
CALL DQRSL (A, LDA, NROWA, NCOLA, QRAUX, Y, NULL, NULL, B,
$ RESID, NULL, JOB, INFO)
IF (INFO .EQ. 0) THEN
PRINT 1040
PRINT 1050, B
PRINT 1060
PRINT 1050, RESID
ELSE
PRINT 1070
END IF
C
1000 FORMAT (1X, 'A:')
1010 FORMAT (2(3X, F4.1))
1020 FORMAT (/1X, 'y:')
1030 FORMAT (3X, F4.1)
1040 FORMAT (/1X, 'Least squares solution:')
1050 FORMAT (3X, F4.1)
1060 FORMAT (/1X, 'Residual:')
1070 FORMAT (1X, 'A is singular.')
C
END
SAMPLE OUTPUT
A:
1.0 1.0
1.0 0.0
0.0 1.0
y:
1.0
0.0
-5.0
Least squares solution:
2.0
-3.0
Residual:
2.0
-2.0
-2.0
|
Закладки на сайте Проследить за страницей |
Created 1996-2025 by Maxim Chirkov Добавить, Поддержать, Вебмастеру |