1      subroutine DGEFSC(A,LDA,N,B,LDB,NB,IPVT,RCOND,Z)
2c Copyright (c) 1996 California Institute of Technology, Pasadena, CA.
3c ALL RIGHTS RESERVED.
4c Based on Government Sponsored Research NAS7-03001.
5C>> 1994-11-11 DGEFSC  Krogh   Declared all vars.
6C>> 1994-10-20 DGEFSC Krogh  Changes to use M77CON
7C>> 1987-08-18 DGEFSC Lawson  Initial code.
8c--D replaces "?": ?GEFSC, ?GECO, ?GESLD
9c
10c     Solves a system of linear equations,  A * X = B,
11c     where A is a square nonsingular matrix of order N and B is an
12c     N by NB matrix.  The solution is the N by NB matrix X that will
13c     be stored on return in place of B in the array B().
14c     Sets RCOND to an estimate of the reciprocal of the condition
15c     number of A.
16c     ------------------------------------------------------------------
17c     Uses subroutines derived from LINPACK.
18c     Ref: LINPACK Users' Guide, by J. J. Dongarra, C. B. Moler,
19c     J. R. Bunch, and G. W. Stewart, publ by Soc. for Indust. and Appl.
20c     Math, Philadelphia, 1979.
21c     Adapted for the JPL Math77 library by C. L. Lawson, JPL, Aug 1987.
22c     ------------------------------------------------------------------
23c                    Subroutine arguments
24c
25c     A(,)  [inout]  On entry contains the N by N matrix A.
26c           On return contains the LU factorization of A as computed by
27c           LINPACK subroutines.
28c
29c     LDA  [in]  Leading dimensioning parameter for the array A(,).
30c           Require LDA .ge. N.
31c
32c     N  [in]  The order of the matrix A and number of rows in the
33c           matrices B and X.
34c
35c     B(,)  [inout]  On entry contains the N by NB matrix, B. The array
36c           B() could be a singly subscripted array if NB .eq. 1.
37c           On return, if RCOND .ne. ZERO, B() will contain the N by NB
38c           solution matrix, X.   If RCOND .eq. ZERO the solution will
39c           not be computed and B() will be unaltered.
40c
41c     LDB  [in]  Leading dimensioning parameter for the array B(,).
42c           Require LDB .ge. N.
43c
44c     NB  [in]  Number of columns in the matrices B and X.  If NB .lt. 1
45c           the matrix, A, will be factored but no reference will be
46c           made to the array B(,).
47c
48c     IPVT()  [out]  Integer array of length at least N.  On return will
49c           contain a record of the row interchanges done during
50c           factorization of A.
51c
52c     RCOND  [out]  On return contains an estimate of the reciprocal of
53c           the condition number of the given matrix, A.  Will be in the
54c           range from zero to one.  Zero indicates an exactly singular
55c           matrix.  A larger RCOND indicates a better conditioned
56c           matrix.  If RCOND .eq. ZERO, the solution X will not be
57c           computed.
58c
59c     Z()  [scratch]  Array of length at least N.  Used as work space.
60c     ------------------------------------------------------------------
61c     Subroutines called: DGECO, DGESLD
62c     ------------------------------------------------------------------
63      integer LDA, N, LDB, NB, IPVT(N), J
64      double precision A(LDA,N), B(LDB,*), Z(N), RCOND, ZERO
65      parameter(ZERO = 0.0D0)
66c     ------------------------------------------------------------------
67c                                    Compute LU factorization of A
68c                                and compute recip. cond. no. of A.
69      call DGECO(A,LDA,N,IPVT,RCOND,Z)
70      if( RCOND .eq. ZERO ) return
71c                                    Solve equations.
72      do 10 J = 1,NB
73         call DGESLD(A,LDA,N,IPVT,B(1,J))
74   10 continue
75      return
76      END
77