Statistiques
| Révision :

root / src / lapack / double / dlascl.f @ 1

Historique | Voir | Annoter | Télécharger (7,76 ko)

1 1 equemene
      SUBROUTINE DLASCL( TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO )
2 1 equemene
*
3 1 equemene
*  -- LAPACK auxiliary routine (version 3.2) --
4 1 equemene
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
5 1 equemene
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
6 1 equemene
*     November 2006
7 1 equemene
*
8 1 equemene
*     .. Scalar Arguments ..
9 1 equemene
      CHARACTER          TYPE
10 1 equemene
      INTEGER            INFO, KL, KU, LDA, M, N
11 1 equemene
      DOUBLE PRECISION   CFROM, CTO
12 1 equemene
*     ..
13 1 equemene
*     .. Array Arguments ..
14 1 equemene
      DOUBLE PRECISION   A( LDA, * )
15 1 equemene
*     ..
16 1 equemene
*
17 1 equemene
*  Purpose
18 1 equemene
*  =======
19 1 equemene
*
20 1 equemene
*  DLASCL multiplies the M by N real matrix A by the real scalar
21 1 equemene
*  CTO/CFROM.  This is done without over/underflow as long as the final
22 1 equemene
*  result CTO*A(I,J)/CFROM does not over/underflow. TYPE specifies that
23 1 equemene
*  A may be full, upper triangular, lower triangular, upper Hessenberg,
24 1 equemene
*  or banded.
25 1 equemene
*
26 1 equemene
*  Arguments
27 1 equemene
*  =========
28 1 equemene
*
29 1 equemene
*  TYPE    (input) CHARACTER*1
30 1 equemene
*          TYPE indices the storage type of the input matrix.
31 1 equemene
*          = 'G':  A is a full matrix.
32 1 equemene
*          = 'L':  A is a lower triangular matrix.
33 1 equemene
*          = 'U':  A is an upper triangular matrix.
34 1 equemene
*          = 'H':  A is an upper Hessenberg matrix.
35 1 equemene
*          = 'B':  A is a symmetric band matrix with lower bandwidth KL
36 1 equemene
*                  and upper bandwidth KU and with the only the lower
37 1 equemene
*                  half stored.
38 1 equemene
*          = 'Q':  A is a symmetric band matrix with lower bandwidth KL
39 1 equemene
*                  and upper bandwidth KU and with the only the upper
40 1 equemene
*                  half stored.
41 1 equemene
*          = 'Z':  A is a band matrix with lower bandwidth KL and upper
42 1 equemene
*                  bandwidth KU.
43 1 equemene
*
44 1 equemene
*  KL      (input) INTEGER
45 1 equemene
*          The lower bandwidth of A.  Referenced only if TYPE = 'B',
46 1 equemene
*          'Q' or 'Z'.
47 1 equemene
*
48 1 equemene
*  KU      (input) INTEGER
49 1 equemene
*          The upper bandwidth of A.  Referenced only if TYPE = 'B',
50 1 equemene
*          'Q' or 'Z'.
51 1 equemene
*
52 1 equemene
*  CFROM   (input) DOUBLE PRECISION
53 1 equemene
*  CTO     (input) DOUBLE PRECISION
54 1 equemene
*          The matrix A is multiplied by CTO/CFROM. A(I,J) is computed
55 1 equemene
*          without over/underflow if the final result CTO*A(I,J)/CFROM
56 1 equemene
*          can be represented without over/underflow.  CFROM must be
57 1 equemene
*          nonzero.
58 1 equemene
*
59 1 equemene
*  M       (input) INTEGER
60 1 equemene
*          The number of rows of the matrix A.  M >= 0.
61 1 equemene
*
62 1 equemene
*  N       (input) INTEGER
63 1 equemene
*          The number of columns of the matrix A.  N >= 0.
64 1 equemene
*
65 1 equemene
*  A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
66 1 equemene
*          The matrix to be multiplied by CTO/CFROM.  See TYPE for the
67 1 equemene
*          storage type.
68 1 equemene
*
69 1 equemene
*  LDA     (input) INTEGER
70 1 equemene
*          The leading dimension of the array A.  LDA >= max(1,M).
71 1 equemene
*
72 1 equemene
*  INFO    (output) INTEGER
73 1 equemene
*          0  - successful exit
74 1 equemene
*          <0 - if INFO = -i, the i-th argument had an illegal value.
75 1 equemene
*
76 1 equemene
*  =====================================================================
77 1 equemene
*
78 1 equemene
*     .. Parameters ..
79 1 equemene
      DOUBLE PRECISION   ZERO, ONE
80 1 equemene
      PARAMETER          ( ZERO = 0.0D0, ONE = 1.0D0 )
81 1 equemene
*     ..
82 1 equemene
*     .. Local Scalars ..
83 1 equemene
      LOGICAL            DONE
84 1 equemene
      INTEGER            I, ITYPE, J, K1, K2, K3, K4
85 1 equemene
      DOUBLE PRECISION   BIGNUM, CFROM1, CFROMC, CTO1, CTOC, MUL, SMLNUM
86 1 equemene
*     ..
87 1 equemene
*     .. External Functions ..
88 1 equemene
      LOGICAL            LSAME, DISNAN
89 1 equemene
      DOUBLE PRECISION   DLAMCH
90 1 equemene
      EXTERNAL           LSAME, DLAMCH, DISNAN
91 1 equemene
*     ..
92 1 equemene
*     .. Intrinsic Functions ..
93 1 equemene
      INTRINSIC          ABS, MAX, MIN
94 1 equemene
*     ..
95 1 equemene
*     .. External Subroutines ..
96 1 equemene
      EXTERNAL           XERBLA
97 1 equemene
*     ..
98 1 equemene
*     .. Executable Statements ..
99 1 equemene
*
100 1 equemene
*     Test the input arguments
101 1 equemene
*
102 1 equemene
      INFO = 0
103 1 equemene
*
104 1 equemene
      IF( LSAME( TYPE, 'G' ) ) THEN
105 1 equemene
         ITYPE = 0
106 1 equemene
      ELSE IF( LSAME( TYPE, 'L' ) ) THEN
107 1 equemene
         ITYPE = 1
108 1 equemene
      ELSE IF( LSAME( TYPE, 'U' ) ) THEN
109 1 equemene
         ITYPE = 2
110 1 equemene
      ELSE IF( LSAME( TYPE, 'H' ) ) THEN
111 1 equemene
         ITYPE = 3
112 1 equemene
      ELSE IF( LSAME( TYPE, 'B' ) ) THEN
113 1 equemene
         ITYPE = 4
114 1 equemene
      ELSE IF( LSAME( TYPE, 'Q' ) ) THEN
115 1 equemene
         ITYPE = 5
116 1 equemene
      ELSE IF( LSAME( TYPE, 'Z' ) ) THEN
117 1 equemene
         ITYPE = 6
118 1 equemene
      ELSE
119 1 equemene
         ITYPE = -1
120 1 equemene
      END IF
121 1 equemene
*
122 1 equemene
      IF( ITYPE.EQ.-1 ) THEN
123 1 equemene
         INFO = -1
124 1 equemene
      ELSE IF( CFROM.EQ.ZERO .OR. DISNAN(CFROM) ) THEN
125 1 equemene
         INFO = -4
126 1 equemene
      ELSE IF( DISNAN(CTO) ) THEN
127 1 equemene
         INFO = -5
128 1 equemene
      ELSE IF( M.LT.0 ) THEN
129 1 equemene
         INFO = -6
130 1 equemene
      ELSE IF( N.LT.0 .OR. ( ITYPE.EQ.4 .AND. N.NE.M ) .OR.
131 1 equemene
     $         ( ITYPE.EQ.5 .AND. N.NE.M ) ) THEN
132 1 equemene
         INFO = -7
133 1 equemene
      ELSE IF( ITYPE.LE.3 .AND. LDA.LT.MAX( 1, M ) ) THEN
134 1 equemene
         INFO = -9
135 1 equemene
      ELSE IF( ITYPE.GE.4 ) THEN
136 1 equemene
         IF( KL.LT.0 .OR. KL.GT.MAX( M-1, 0 ) ) THEN
137 1 equemene
            INFO = -2
138 1 equemene
         ELSE IF( KU.LT.0 .OR. KU.GT.MAX( N-1, 0 ) .OR.
139 1 equemene
     $            ( ( ITYPE.EQ.4 .OR. ITYPE.EQ.5 ) .AND. KL.NE.KU ) )
140 1 equemene
     $             THEN
141 1 equemene
            INFO = -3
142 1 equemene
         ELSE IF( ( ITYPE.EQ.4 .AND. LDA.LT.KL+1 ) .OR.
143 1 equemene
     $            ( ITYPE.EQ.5 .AND. LDA.LT.KU+1 ) .OR.
144 1 equemene
     $            ( ITYPE.EQ.6 .AND. LDA.LT.2*KL+KU+1 ) ) THEN
145 1 equemene
            INFO = -9
146 1 equemene
         END IF
147 1 equemene
      END IF
148 1 equemene
*
149 1 equemene
      IF( INFO.NE.0 ) THEN
150 1 equemene
         CALL XERBLA( 'DLASCL', -INFO )
151 1 equemene
         RETURN
152 1 equemene
      END IF
153 1 equemene
*
154 1 equemene
*     Quick return if possible
155 1 equemene
*
156 1 equemene
      IF( N.EQ.0 .OR. M.EQ.0 )
157 1 equemene
     $   RETURN
158 1 equemene
*
159 1 equemene
*     Get machine parameters
160 1 equemene
*
161 1 equemene
      SMLNUM = DLAMCH( 'S' )
162 1 equemene
      BIGNUM = ONE / SMLNUM
163 1 equemene
*
164 1 equemene
      CFROMC = CFROM
165 1 equemene
      CTOC = CTO
166 1 equemene
*
167 1 equemene
   10 CONTINUE
168 1 equemene
      CFROM1 = CFROMC*SMLNUM
169 1 equemene
      IF( CFROM1.EQ.CFROMC ) THEN
170 1 equemene
!        CFROMC is an inf.  Multiply by a correctly signed zero for
171 1 equemene
!        finite CTOC, or a NaN if CTOC is infinite.
172 1 equemene
         MUL = CTOC / CFROMC
173 1 equemene
         DONE = .TRUE.
174 1 equemene
         CTO1 = CTOC
175 1 equemene
      ELSE
176 1 equemene
         CTO1 = CTOC / BIGNUM
177 1 equemene
         IF( CTO1.EQ.CTOC ) THEN
178 1 equemene
!           CTOC is either 0 or an inf.  In both cases, CTOC itself
179 1 equemene
!           serves as the correct multiplication factor.
180 1 equemene
            MUL = CTOC
181 1 equemene
            DONE = .TRUE.
182 1 equemene
            CFROMC = ONE
183 1 equemene
         ELSE IF( ABS( CFROM1 ).GT.ABS( CTOC ) .AND. CTOC.NE.ZERO ) THEN
184 1 equemene
            MUL = SMLNUM
185 1 equemene
            DONE = .FALSE.
186 1 equemene
            CFROMC = CFROM1
187 1 equemene
         ELSE IF( ABS( CTO1 ).GT.ABS( CFROMC ) ) THEN
188 1 equemene
            MUL = BIGNUM
189 1 equemene
            DONE = .FALSE.
190 1 equemene
            CTOC = CTO1
191 1 equemene
         ELSE
192 1 equemene
            MUL = CTOC / CFROMC
193 1 equemene
            DONE = .TRUE.
194 1 equemene
         END IF
195 1 equemene
      END IF
196 1 equemene
*
197 1 equemene
      IF( ITYPE.EQ.0 ) THEN
198 1 equemene
*
199 1 equemene
*        Full matrix
200 1 equemene
*
201 1 equemene
         DO 30 J = 1, N
202 1 equemene
            DO 20 I = 1, M
203 1 equemene
               A( I, J ) = A( I, J )*MUL
204 1 equemene
   20       CONTINUE
205 1 equemene
   30    CONTINUE
206 1 equemene
*
207 1 equemene
      ELSE IF( ITYPE.EQ.1 ) THEN
208 1 equemene
*
209 1 equemene
*        Lower triangular matrix
210 1 equemene
*
211 1 equemene
         DO 50 J = 1, N
212 1 equemene
            DO 40 I = J, M
213 1 equemene
               A( I, J ) = A( I, J )*MUL
214 1 equemene
   40       CONTINUE
215 1 equemene
   50    CONTINUE
216 1 equemene
*
217 1 equemene
      ELSE IF( ITYPE.EQ.2 ) THEN
218 1 equemene
*
219 1 equemene
*        Upper triangular matrix
220 1 equemene
*
221 1 equemene
         DO 70 J = 1, N
222 1 equemene
            DO 60 I = 1, MIN( J, M )
223 1 equemene
               A( I, J ) = A( I, J )*MUL
224 1 equemene
   60       CONTINUE
225 1 equemene
   70    CONTINUE
226 1 equemene
*
227 1 equemene
      ELSE IF( ITYPE.EQ.3 ) THEN
228 1 equemene
*
229 1 equemene
*        Upper Hessenberg matrix
230 1 equemene
*
231 1 equemene
         DO 90 J = 1, N
232 1 equemene
            DO 80 I = 1, MIN( J+1, M )
233 1 equemene
               A( I, J ) = A( I, J )*MUL
234 1 equemene
   80       CONTINUE
235 1 equemene
   90    CONTINUE
236 1 equemene
*
237 1 equemene
      ELSE IF( ITYPE.EQ.4 ) THEN
238 1 equemene
*
239 1 equemene
*        Lower half of a symmetric band matrix
240 1 equemene
*
241 1 equemene
         K3 = KL + 1
242 1 equemene
         K4 = N + 1
243 1 equemene
         DO 110 J = 1, N
244 1 equemene
            DO 100 I = 1, MIN( K3, K4-J )
245 1 equemene
               A( I, J ) = A( I, J )*MUL
246 1 equemene
  100       CONTINUE
247 1 equemene
  110    CONTINUE
248 1 equemene
*
249 1 equemene
      ELSE IF( ITYPE.EQ.5 ) THEN
250 1 equemene
*
251 1 equemene
*        Upper half of a symmetric band matrix
252 1 equemene
*
253 1 equemene
         K1 = KU + 2
254 1 equemene
         K3 = KU + 1
255 1 equemene
         DO 130 J = 1, N
256 1 equemene
            DO 120 I = MAX( K1-J, 1 ), K3
257 1 equemene
               A( I, J ) = A( I, J )*MUL
258 1 equemene
  120       CONTINUE
259 1 equemene
  130    CONTINUE
260 1 equemene
*
261 1 equemene
      ELSE IF( ITYPE.EQ.6 ) THEN
262 1 equemene
*
263 1 equemene
*        Band matrix
264 1 equemene
*
265 1 equemene
         K1 = KL + KU + 2
266 1 equemene
         K2 = KL + 1
267 1 equemene
         K3 = 2*KL + KU + 1
268 1 equemene
         K4 = KL + KU + 1 + M
269 1 equemene
         DO 150 J = 1, N
270 1 equemene
            DO 140 I = MAX( K1-J, K2 ), MIN( K3, K4-J )
271 1 equemene
               A( I, J ) = A( I, J )*MUL
272 1 equemene
  140       CONTINUE
273 1 equemene
  150    CONTINUE
274 1 equemene
*
275 1 equemene
      END IF
276 1 equemene
*
277 1 equemene
      IF( .NOT.DONE )
278 1 equemene
     $   GO TO 10
279 1 equemene
*
280 1 equemene
      RETURN
281 1 equemene
*
282 1 equemene
*     End of DLASCL
283 1 equemene
*
284 1 equemene
      END