Télécharger dtrsyl.eso

Retour à la liste

Numérotation des lignes :

dtrsyl
  1. C DTRSYL SOURCE FANDEUR 22/05/02 21:15:18 11359
  2. *> \brief \b DTRSYL
  3. *
  4. * =========== DOCUMENTATION ===========
  5. *
  6. * Online html documentation available at
  7. * http://www.netlib.org/lapack/explore-html/
  8. *
  9. *> \htmlonly
  10. *> Download DTRSYL + dependencies
  11. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dtrsyl.f">
  12. *> [TGZ]</a>
  13. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dtrsyl.f">
  14. *> [ZIP]</a>
  15. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dtrsyl.f">
  16. *> [TXT]</a>
  17. *> \endhtmlonly
  18. *
  19. * Definition:
  20. * ===========
  21. *
  22. * SUBROUTINE DTRSYL( TRANA, TRANB, ISGN, M, N, A, LDA, B, LDB, C,
  23. * LDC, SCALE, INFO )
  24. *
  25. * .. Scalar Arguments ..
  26. * CHARACTER TRANA, TRANB
  27. * INTEGER INFO, ISGN, LDA, LDB, LDC, M, N
  28. * REAL*8 SCALE
  29. * ..
  30. * .. Array Arguments ..
  31. * REAL*8 A( LDA, * ), B( LDB, * ), C( LDC, * )
  32. * ..
  33. *
  34. *
  35. *> \par Purpose:
  36. * =============
  37. *>
  38. *> \verbatim
  39. *>
  40. *> DTRSYL solves the real Sylvester matrix equation:
  41. *>
  42. *> op(A)*X + X*op(B) = scale*C or
  43. *> op(A)*X - X*op(B) = scale*C,
  44. *>
  45. *> where op(A) = A or A**T, and A and B are both upper quasi-
  46. *> triangular. A is M-by-M and B is N-by-N; the right hand side C and
  47. *> the solution X are M-by-N; and scale is an output scale factor, set
  48. *> <= 1 to avoid overflow in X.
  49. *>
  50. *> A and B must be in Schur canonical form (as returned by DHSEQR), that
  51. *> is, block upper triangular with 1-by-1 and 2-by-2 diagonal blocks;
  52. *> each 2-by-2 diagonal block has its diagonal elements equal and its
  53. *> off-diagonal elements of opposite sign.
  54. *> \endverbatim
  55. *
  56. * Arguments:
  57. * ==========
  58. *
  59. *> \param[in] TRANA
  60. *> \verbatim
  61. *> TRANA is CHARACTER*1
  62. *> Specifies the option op(A):
  63. *> = 'N': op(A) = A (No transpose)
  64. *> = 'T': op(A) = A**T (Transpose)
  65. *> = 'C': op(A) = A**H (Conjugate transpose = Transpose)
  66. *> \endverbatim
  67. *>
  68. *> \param[in] TRANB
  69. *> \verbatim
  70. *> TRANB is CHARACTER*1
  71. *> Specifies the option op(B):
  72. *> = 'N': op(B) = B (No transpose)
  73. *> = 'T': op(B) = B**T (Transpose)
  74. *> = 'C': op(B) = B**H (Conjugate transpose = Transpose)
  75. *> \endverbatim
  76. *>
  77. *> \param[in] ISGN
  78. *> \verbatim
  79. *> ISGN is INTEGER
  80. *> Specifies the sign in the equation:
  81. *> = +1: solve op(A)*X + X*op(B) = scale*C
  82. *> = -1: solve op(A)*X - X*op(B) = scale*C
  83. *> \endverbatim
  84. *>
  85. *> \param[in] M
  86. *> \verbatim
  87. *> M is INTEGER
  88. *> The order of the matrix A, and the number of rows in the
  89. *> matrices X and C. M >= 0.
  90. *> \endverbatim
  91. *>
  92. *> \param[in] N
  93. *> \verbatim
  94. *> N is INTEGER
  95. *> The order of the matrix B, and the number of columns in the
  96. *> matrices X and C. N >= 0.
  97. *> \endverbatim
  98. *>
  99. *> \param[in] A
  100. *> \verbatim
  101. *> A is DOUBLE PRECISION array, dimension (LDA,M)
  102. *> The upper quasi-triangular matrix A, in Schur canonical form.
  103. *> \endverbatim
  104. *>
  105. *> \param[in] LDA
  106. *> \verbatim
  107. *> LDA is INTEGER
  108. *> The leading dimension of the array A. LDA >= max(1,M).
  109. *> \endverbatim
  110. *>
  111. *> \param[in] B
  112. *> \verbatim
  113. *> B is DOUBLE PRECISION array, dimension (LDB,N)
  114. *> The upper quasi-triangular matrix B, in Schur canonical form.
  115. *> \endverbatim
  116. *>
  117. *> \param[in] LDB
  118. *> \verbatim
  119. *> LDB is INTEGER
  120. *> The leading dimension of the array B. LDB >= max(1,N).
  121. *> \endverbatim
  122. *>
  123. *> \param[in,out] C
  124. *> \verbatim
  125. *> C is DOUBLE PRECISION array, dimension (LDC,N)
  126. *> On entry, the M-by-N right hand side matrix C.
  127. *> On exit, C is overwritten by the solution matrix X.
  128. *> \endverbatim
  129. *>
  130. *> \param[in] LDC
  131. *> \verbatim
  132. *> LDC is INTEGER
  133. *> The leading dimension of the array C. LDC >= max(1,M)
  134. *> \endverbatim
  135. *>
  136. *> \param[out] SCALE
  137. *> \verbatim
  138. *> SCALE is DOUBLE PRECISION
  139. *> The scale factor, scale, set <= 1 to avoid overflow in X.
  140. *> \endverbatim
  141. *>
  142. *> \param[out] INFO
  143. *> \verbatim
  144. *> INFO is INTEGER
  145. *> = 0: successful exit
  146. *> < 0: if INFO = -i, the i-th argument had an illegal value
  147. *> = 1: A and B have common or very close eigenvalues; perturbed
  148. *> values were used to solve the equation (but the matrices
  149. *> A and B are unchanged).
  150. *> \endverbatim
  151. *
  152. * Authors:
  153. * ========
  154. *
  155. *> \author Univ. of Tennessee
  156. *> \author Univ. of California Berkeley
  157. *> \author Univ. of Colorado Denver
  158. *> \author NAG Ltd.
  159. *
  160. *> \date December 2016
  161. *
  162. *> \ingroup doubleSYcomputational
  163. *
  164. * =====================================================================
  165. SUBROUTINE DTRSYL( TRANA, TRANB, ISGN, M, N, A, LDA, B, LDB, C,
  166. $ LDC, SCALE, INFO )
  167. *
  168. * -- LAPACK computational routine (version 3.7.0) --
  169. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  170. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  171. * December 2016
  172. *
  173. * .. Scalar Arguments ..
  174. CHARACTER TRANA, TRANB
  175. INTEGER INFO, ISGN, LDA, LDB, LDC, M, N
  176. REAL*8 SCALE
  177. * ..
  178. * .. Array Arguments ..
  179. REAL*8 A( LDA, * ), B( LDB, * ), C( LDC, * )
  180. * ..
  181. *
  182. * =====================================================================
  183. *
  184. * .. Parameters ..
  185. REAL*8 ZERO, ONE
  186. PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
  187. * ..
  188. * .. Local Scalars ..
  189. LOGICAL NOTRNA, NOTRNB
  190. INTEGER IERR, J, K, K1, K2, KNEXT, L, L1, L2, LNEXT
  191. REAL*8 A11, BIGNUM, DA11, DB, EPS, SCALOC, SGN, SMIN,
  192. $ SMLNUM, SUML, SUMR, XNORM
  193. * ..
  194. * .. Local Arrays ..
  195. REAL*8 DUM( 1 ), VEC( 2, 2 ), X( 2, 2 )
  196. * ..
  197. * .. External Functions ..
  198. LOGICAL LSAME
  199. REAL*8 DDOT, DLAMCH, DLANGE
  200. EXTERNAL LSAME, DDOT, DLAMCH, DLANGE
  201. * ..
  202. * .. External Subroutines ..
  203. EXTERNAL DLABAD, DLALN2, DLASY2,
  204. * ..
  205. ** .. Intrinsic Functions ..
  206. * INTRINSIC ABS, DBLE, MAX, MIN
  207. ** ..
  208. ** .. Executable Statements ..
  209. *
  210. * Decode and Test input parameters
  211. *
  212. NOTRNA = LSAME( TRANA, 'N' )
  213. NOTRNB = LSAME( TRANB, 'N' )
  214. *
  215. INFO = 0
  216. IF( .NOT.NOTRNA .AND.
  217. $ .NOT.LSAME( TRANA, 'T' ) .AND.
  218. $ .NOT.LSAME( TRANA, 'C' ) ) THEN
  219. INFO = -1
  220. ELSE IF( .NOT.NOTRNB .AND.
  221. $ .NOT.LSAME( TRANB, 'T' ) .AND.
  222. $ .NOT.LSAME( TRANB, 'C' ) ) THEN
  223. INFO = -2
  224. ELSE IF( ISGN.NE.1 .AND. ISGN.NE.-1 ) THEN
  225. INFO = -3
  226. ELSE IF( M.LT.0 ) THEN
  227. INFO = -4
  228. ELSE IF( N.LT.0 ) THEN
  229. INFO = -5
  230. ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
  231. INFO = -7
  232. ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
  233. INFO = -9
  234. ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
  235. INFO = -11
  236. END IF
  237. IF( INFO.NE.0 ) THEN
  238. CALL XERBLA( 'DTRSYL', -INFO )
  239. RETURN
  240. END IF
  241. *
  242. * Quick return if possible
  243. *
  244. SCALE = ONE
  245. IF( M.EQ.0 .OR. N.EQ.0 )
  246. $ RETURN
  247. *
  248. * Set constants to control overflow
  249. *
  250. EPS = DLAMCH( 'P' )
  251. SMLNUM = DLAMCH( 'S' )
  252. BIGNUM = ONE / SMLNUM
  253. CALL DLABAD( SMLNUM, BIGNUM )
  254. SMLNUM = SMLNUM*DBLE( M*N ) / EPS
  255. BIGNUM = ONE / SMLNUM
  256. *
  257. SMIN = MAX( SMLNUM, EPS*DLANGE( 'M', M, M, A, LDA, DUM ),
  258. $ EPS*DLANGE( 'M', N, N, B, LDB, DUM ) )
  259. *
  260. SGN = ISGN
  261. *
  262. IF( NOTRNA .AND. NOTRNB ) THEN
  263. *
  264. * Solve A*X + ISGN*X*B = scale*C.
  265. *
  266. * The (K,L)th block of X is determined starting from
  267. * bottom-left corner column by column by
  268. *
  269. * A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
  270. *
  271. * Where
  272. * M L-1
  273. * R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(J,L)].
  274. * I=K+1 J=1
  275. *
  276. * Start column loop (index = L)
  277. * L1 (L2) : column index of the first (first) row of X(K,L).
  278. *
  279. LNEXT = 1
  280. DO 60 L = 1, N
  281. IF( L.LT.LNEXT )
  282. $ GO TO 60
  283. IF( L.EQ.N ) THEN
  284. L1 = L
  285. L2 = L
  286. ELSE
  287. IF( B( L+1, L ).NE.ZERO ) THEN
  288. L1 = L
  289. L2 = L + 1
  290. LNEXT = L + 2
  291. ELSE
  292. L1 = L
  293. L2 = L
  294. LNEXT = L + 1
  295. END IF
  296. END IF
  297. *
  298. * Start row loop (index = K)
  299. * K1 (K2): row index of the first (last) row of X(K,L).
  300. *
  301. KNEXT = M
  302. DO 50 K = M, 1, -1
  303. IF( K.GT.KNEXT )
  304. $ GO TO 50
  305. IF( K.EQ.1 ) THEN
  306. K1 = K
  307. K2 = K
  308. ELSE
  309. IF( A( K, K-1 ).NE.ZERO ) THEN
  310. K1 = K - 1
  311. K2 = K
  312. KNEXT = K - 2
  313. ELSE
  314. K1 = K
  315. K2 = K
  316. KNEXT = K - 1
  317. END IF
  318. END IF
  319. *
  320. IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
  321. SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ),
  322. $ LDA, C( MIN( K1+1, M ), L1 ), 1 )
  323. SUMR = DDOT( L1-1, C( K1, 1 ),
  324. $ LDC, B( 1, L1 ), 1 )
  325. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  326. SCALOC = ONE
  327. *
  328. A11 = A( K1, K1 ) + SGN*B( L1, L1 )
  329. DA11 = ABS( A11 )
  330. IF( DA11.LE.SMIN ) THEN
  331. A11 = SMIN
  332. DA11 = SMIN
  333. INFO = 1
  334. END IF
  335. DB = ABS( VEC( 1, 1 ) )
  336. IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
  337. IF( DB.GT.BIGNUM*DA11 )
  338. $ SCALOC = ONE / DB
  339. END IF
  340. X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
  341. *
  342. IF( SCALOC.NE.ONE ) THEN
  343. DO 10 J = 1, N
  344. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  345. 10 CONTINUE
  346. SCALE = SCALE*SCALOC
  347. END IF
  348. C( K1, L1 ) = X( 1, 1 )
  349. *
  350. ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
  351. *
  352. SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ),
  353. $ LDA, C( MIN( K2+1, M ), L1 ), 1 )
  354. SUMR = DDOT( L1-1, C( K1, 1 ),
  355. $ LDC, B( 1, L1 ), 1 )
  356. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  357. *
  358. SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ),
  359. $ LDA, C( MIN( K2+1, M ), L1 ), 1 )
  360. SUMR = DDOT( L1-1, C( K2, 1 ),
  361. $ LDC, B( 1, L1 ), 1 )
  362. VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
  363. *
  364. CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, A( K1, K1 ),
  365. $ LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
  366. $ ZERO, X, 2, SCALOC, XNORM, IERR )
  367. IF( IERR.NE.0 )
  368. $ INFO = 1
  369. *
  370. IF( SCALOC.NE.ONE ) THEN
  371. DO 20 J = 1, N
  372. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  373. 20 CONTINUE
  374. SCALE = SCALE*SCALOC
  375. END IF
  376. C( K1, L1 ) = X( 1, 1 )
  377. C( K2, L1 ) = X( 2, 1 )
  378. *
  379. ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
  380. *
  381. SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ),
  382. $ LDA, C( MIN( K1+1, M ), L1 ), 1 )
  383. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  384. $ B( 1, L1 ), 1 )
  385. VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
  386. *
  387. SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ),
  388. $ LDA, C( MIN( K1+1, M ), L2 ), 1 )
  389. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  390. $ B( 1, L2 ), 1 )
  391. VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
  392. *
  393. CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, B( L1, L1 ),
  394. $ LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
  395. $ ZERO, X, 2, SCALOC, XNORM, IERR )
  396. IF( IERR.NE.0 )
  397. $ INFO = 1
  398. *
  399. IF( SCALOC.NE.ONE ) THEN
  400. DO 30 J = 1, N
  401. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  402. 30 CONTINUE
  403. SCALE = SCALE*SCALOC
  404. END IF
  405. C( K1, L1 ) = X( 1, 1 )
  406. C( K1, L2 ) = X( 2, 1 )
  407. *
  408. ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
  409. *
  410. SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ),
  411. $ LDA, C( MIN( K2+1, M ), L1 ), 1 )
  412. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  413. $ B( 1, L1 ), 1 )
  414. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  415. *
  416. SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ),
  417. $ LDA, C( MIN( K2+1, M ), L2 ), 1 )
  418. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  419. $ B( 1, L2 ), 1 )
  420. VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
  421. *
  422. SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ),
  423. $ LDA, C( MIN( K2+1, M ), L1 ), 1 )
  424. SUMR = DDOT( L1-1, C( K2, 1 ), LDC,
  425. $ B( 1, L1 ), 1 )
  426. VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
  427. *
  428. SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ),
  429. $ LDA, C( MIN( K2+1, M ), L2 ), 1 )
  430. SUMR = DDOT( L1-1, C( K2, 1 ), LDC,
  431. $ B( 1, L2 ), 1 )
  432. VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
  433. *
  434. CALL DLASY2( .FALSE., .FALSE., ISGN, 2, 2,
  435. $ A( K1, K1 ), LDA, B( L1, L1 ), LDB, VEC,
  436. $ 2, SCALOC, X, 2, XNORM, IERR )
  437. IF( IERR.NE.0 )
  438. $ INFO = 1
  439. *
  440. IF( SCALOC.NE.ONE ) THEN
  441. DO 40 J = 1, N
  442. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  443. 40 CONTINUE
  444. SCALE = SCALE*SCALOC
  445. END IF
  446. C( K1, L1 ) = X( 1, 1 )
  447. C( K1, L2 ) = X( 1, 2 )
  448. C( K2, L1 ) = X( 2, 1 )
  449. C( K2, L2 ) = X( 2, 2 )
  450. END IF
  451. *
  452. 50 CONTINUE
  453. *
  454. 60 CONTINUE
  455. *
  456. ELSE IF( .NOT.NOTRNA .AND. NOTRNB ) THEN
  457. *
  458. * Solve A**T *X + ISGN*X*B = scale*C.
  459. *
  460. * The (K,L)th block of X is determined starting from
  461. * upper-left corner column by column by
  462. *
  463. * A(K,K)**T*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
  464. *
  465. * Where
  466. * K-1 T L-1
  467. * R(K,L) = SUM [A(I,K)**T*X(I,L)] +ISGN*SUM [X(K,J)*B(J,L)]
  468. * I=1 J=1
  469. *
  470. * Start column loop (index = L)
  471. * L1 (L2): column index of the first (last) row of X(K,L)
  472. *
  473. LNEXT = 1
  474. DO 120 L = 1, N
  475. IF( L.LT.LNEXT )
  476. $ GO TO 120
  477. IF( L.EQ.N ) THEN
  478. L1 = L
  479. L2 = L
  480. ELSE
  481. IF( B( L+1, L ).NE.ZERO ) THEN
  482. L1 = L
  483. L2 = L + 1
  484. LNEXT = L + 2
  485. ELSE
  486. L1 = L
  487. L2 = L
  488. LNEXT = L + 1
  489. END IF
  490. END IF
  491. *
  492. * Start row loop (index = K)
  493. * K1 (K2): row index of the first (last) row of X(K,L)
  494. *
  495. KNEXT = 1
  496. DO 110 K = 1, M
  497. IF( K.LT.KNEXT )
  498. $ GO TO 110
  499. IF( K.EQ.M ) THEN
  500. K1 = K
  501. K2 = K
  502. ELSE
  503. IF( A( K+1, K ).NE.ZERO ) THEN
  504. K1 = K
  505. K2 = K + 1
  506. KNEXT = K + 2
  507. ELSE
  508. K1 = K
  509. K2 = K
  510. KNEXT = K + 1
  511. END IF
  512. END IF
  513. *
  514. IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
  515. SUML = DDOT( K1-1, A( 1, K1 ), 1,
  516. $ C( 1, L1 ), 1 )
  517. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  518. $ B( 1, L1 ), 1 )
  519. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  520. SCALOC = ONE
  521. *
  522. A11 = A( K1, K1 ) + SGN*B( L1, L1 )
  523. DA11 = ABS( A11 )
  524. IF( DA11.LE.SMIN ) THEN
  525. A11 = SMIN
  526. DA11 = SMIN
  527. INFO = 1
  528. END IF
  529. DB = ABS( VEC( 1, 1 ) )
  530. IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
  531. IF( DB.GT.BIGNUM*DA11 )
  532. $ SCALOC = ONE / DB
  533. END IF
  534. X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
  535. *
  536. IF( SCALOC.NE.ONE ) THEN
  537. DO 70 J = 1, N
  538. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  539. 70 CONTINUE
  540. SCALE = SCALE*SCALOC
  541. END IF
  542. C( K1, L1 ) = X( 1, 1 )
  543. *
  544. ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
  545. *
  546. SUML = DDOT( K1-1, A( 1, K1 ), 1,
  547. $ C( 1, L1 ), 1 )
  548. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  549. $ B( 1, L1 ), 1 )
  550. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  551. *
  552. SUML = DDOT( K1-1, A( 1, K2 ), 1,
  553. $ C( 1, L1 ), 1 )
  554. SUMR = DDOT( L1-1, C( K2, 1 ), LDC,
  555. $ B( 1, L1 ), 1 )
  556. VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
  557. *
  558. CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, A( K1, K1 ),
  559. $ LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
  560. $ ZERO, X, 2, SCALOC, XNORM, IERR )
  561. IF( IERR.NE.0 )
  562. $ INFO = 1
  563. *
  564. IF( SCALOC.NE.ONE ) THEN
  565. DO 80 J = 1, N
  566. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  567. 80 CONTINUE
  568. SCALE = SCALE*SCALOC
  569. END IF
  570. C( K1, L1 ) = X( 1, 1 )
  571. C( K2, L1 ) = X( 2, 1 )
  572. *
  573. ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
  574. *
  575. SUML = DDOT( K1-1, A( 1, K1 ), 1,
  576. $ C( 1, L1 ), 1 )
  577. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  578. $ B( 1, L1 ), 1 )
  579. VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
  580. *
  581. SUML = DDOT( K1-1, A( 1, K1 ), 1,
  582. $ C( 1, L2 ), 1 )
  583. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  584. $ B( 1, L2 ), 1 )
  585. VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
  586. *
  587. CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, B( L1, L1 ),
  588. $ LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
  589. $ ZERO, X, 2, SCALOC, XNORM, IERR )
  590. IF( IERR.NE.0 )
  591. $ INFO = 1
  592. *
  593. IF( SCALOC.NE.ONE ) THEN
  594. DO 90 J = 1, N
  595. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  596. 90 CONTINUE
  597. SCALE = SCALE*SCALOC
  598. END IF
  599. C( K1, L1 ) = X( 1, 1 )
  600. C( K1, L2 ) = X( 2, 1 )
  601. *
  602. ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
  603. *
  604. SUML = DDOT( K1-1, A( 1, K1 ), 1,
  605. $ C( 1, L1 ), 1 )
  606. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  607. $ B( 1, L1 ), 1 )
  608. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  609. *
  610. SUML = DDOT( K1-1, A( 1, K1 ), 1,
  611. $ C( 1, L2 ), 1 )
  612. SUMR = DDOT( L1-1, C( K1, 1 ), LDC,
  613. $ B( 1, L2 ), 1 )
  614. VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
  615. *
  616. SUML = DDOT( K1-1, A( 1, K2 ), 1,
  617. $ C( 1, L1 ), 1 )
  618. SUMR = DDOT( L1-1, C( K2, 1 ), LDC,
  619. $ B( 1, L1 ), 1 )
  620. VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
  621. *
  622. SUML = DDOT( K1-1, A( 1, K2 ), 1,
  623. $ C( 1, L2 ), 1 )
  624. SUMR = DDOT( L1-1, C( K2, 1 ), LDC,
  625. $ B( 1, L2 ), 1 )
  626. VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
  627. *
  628. CALL DLASY2( .TRUE., .FALSE., ISGN, 2, 2, A( K1, K1 ),
  629. $ LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
  630. $ 2, XNORM, IERR )
  631. IF( IERR.NE.0 )
  632. $ INFO = 1
  633. *
  634. IF( SCALOC.NE.ONE ) THEN
  635. DO 100 J = 1, N
  636. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  637. 100 CONTINUE
  638. SCALE = SCALE*SCALOC
  639. END IF
  640. C( K1, L1 ) = X( 1, 1 )
  641. C( K1, L2 ) = X( 1, 2 )
  642. C( K2, L1 ) = X( 2, 1 )
  643. C( K2, L2 ) = X( 2, 2 )
  644. END IF
  645. *
  646. 110 CONTINUE
  647. 120 CONTINUE
  648. *
  649. ELSE IF( .NOT.NOTRNA .AND. .NOT.NOTRNB ) THEN
  650. *
  651. * Solve A**T*X + ISGN*X*B**T = scale*C.
  652. *
  653. * The (K,L)th block of X is determined starting from
  654. * top-right corner column by column by
  655. *
  656. * A(K,K)**T*X(K,L) + ISGN*X(K,L)*B(L,L)**T = C(K,L) - R(K,L)
  657. *
  658. * Where
  659. * K-1 N
  660. * R(K,L) = SUM [A(I,K)**T*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)**T].
  661. * I=1 J=L+1
  662. *
  663. * Start column loop (index = L)
  664. * L1 (L2): column index of the first (last) row of X(K,L)
  665. *
  666. LNEXT = N
  667. DO 180 L = N, 1, -1
  668. IF( L.GT.LNEXT )
  669. $ GO TO 180
  670. IF( L.EQ.1 ) THEN
  671. L1 = L
  672. L2 = L
  673. ELSE
  674. IF( B( L, L-1 ).NE.ZERO ) THEN
  675. L1 = L - 1
  676. L2 = L
  677. LNEXT = L - 2
  678. ELSE
  679. L1 = L
  680. L2 = L
  681. LNEXT = L - 1
  682. END IF
  683. END IF
  684. *
  685. * Start row loop (index = K)
  686. * K1 (K2): row index of the first (last) row of X(K,L)
  687. *
  688. KNEXT = 1
  689. DO 170 K = 1, M
  690. IF( K.LT.KNEXT )
  691. $ GO TO 170
  692. IF( K.EQ.M ) THEN
  693. K1 = K
  694. K2 = K
  695. ELSE
  696. IF( A( K+1, K ).NE.ZERO ) THEN
  697. K1 = K
  698. K2 = K + 1
  699. KNEXT = K + 2
  700. ELSE
  701. K1 = K
  702. K2 = K
  703. KNEXT = K + 1
  704. END IF
  705. END IF
  706. *
  707. IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
  708. SUML = DDOT( K1-1, A( 1, K1 ), 1,
  709. $ C( 1, L1 ), 1 )
  710. SUMR = DDOT( N-L1, C( K1, MIN( L1+1, N ) ),
  711. $ LDC, B( L1, MIN( L1+1, N ) ), LDB )
  712. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  713. SCALOC = ONE
  714. *
  715. A11 = A( K1, K1 ) + SGN*B( L1, L1 )
  716. DA11 = ABS( A11 )
  717. IF( DA11.LE.SMIN ) THEN
  718. A11 = SMIN
  719. DA11 = SMIN
  720. INFO = 1
  721. END IF
  722. DB = ABS( VEC( 1, 1 ) )
  723. IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
  724. IF( DB.GT.BIGNUM*DA11 )
  725. $ SCALOC = ONE / DB
  726. END IF
  727. X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
  728. *
  729. IF( SCALOC.NE.ONE ) THEN
  730. DO 130 J = 1, N
  731. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  732. 130 CONTINUE
  733. SCALE = SCALE*SCALOC
  734. END IF
  735. C( K1, L1 ) = X( 1, 1 )
  736. *
  737. ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
  738. *
  739. SUML = DDOT( K1-1, A( 1, K1 ), 1,
  740. $ C( 1, L1 ), 1 )
  741. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  742. $ B( L1, MIN( L2+1, N ) ), LDB )
  743. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  744. *
  745. SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
  746. SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
  747. $ B( L1, MIN( L2+1, N ) ), LDB )
  748. VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
  749. *
  750. CALL DLALN2( .TRUE., 2, 1, SMIN, ONE, A( K1, K1 ),
  751. $ LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
  752. $ ZERO, X, 2, SCALOC, XNORM, IERR )
  753. IF( IERR.NE.0 )
  754. $ INFO = 1
  755. *
  756. IF( SCALOC.NE.ONE ) THEN
  757. DO 140 J = 1, N
  758. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  759. 140 CONTINUE
  760. SCALE = SCALE*SCALOC
  761. END IF
  762. C( K1, L1 ) = X( 1, 1 )
  763. C( K2, L1 ) = X( 2, 1 )
  764. *
  765. ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
  766. *
  767. SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
  768. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  769. $ B( L1, MIN( L2+1, N ) ), LDB )
  770. VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
  771. *
  772. SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
  773. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  774. $ B( L2, MIN( L2+1, N ) ), LDB )
  775. VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
  776. *
  777. CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, B( L1, L1 ),
  778. $ LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
  779. $ ZERO, X, 2, SCALOC, XNORM, IERR )
  780. IF( IERR.NE.0 )
  781. $ INFO = 1
  782. *
  783. IF( SCALOC.NE.ONE ) THEN
  784. DO 150 J = 1, N
  785. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  786. 150 CONTINUE
  787. SCALE = SCALE*SCALOC
  788. END IF
  789. C( K1, L1 ) = X( 1, 1 )
  790. C( K1, L2 ) = X( 2, 1 )
  791. *
  792. ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
  793. *
  794. SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 )
  795. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  796. $ B( L1, MIN( L2+1, N ) ), LDB )
  797. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  798. *
  799. SUML = DDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 )
  800. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  801. $ B( L2, MIN( L2+1, N ) ), LDB )
  802. VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
  803. *
  804. SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 )
  805. SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
  806. $ B( L1, MIN( L2+1, N ) ), LDB )
  807. VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
  808. *
  809. SUML = DDOT( K1-1, A( 1, K2 ), 1, C( 1, L2 ), 1 )
  810. SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
  811. $ B( L2, MIN( L2+1, N ) ), LDB )
  812. VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
  813. *
  814. CALL DLASY2( .TRUE., .TRUE., ISGN, 2, 2, A( K1, K1 ),
  815. $ LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
  816. $ 2, XNORM, IERR )
  817. IF( IERR.NE.0 )
  818. $ INFO = 1
  819. *
  820. IF( SCALOC.NE.ONE ) THEN
  821. DO 160 J = 1, N
  822. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  823. 160 CONTINUE
  824. SCALE = SCALE*SCALOC
  825. END IF
  826. C( K1, L1 ) = X( 1, 1 )
  827. C( K1, L2 ) = X( 1, 2 )
  828. C( K2, L1 ) = X( 2, 1 )
  829. C( K2, L2 ) = X( 2, 2 )
  830. END IF
  831. *
  832. 170 CONTINUE
  833. 180 CONTINUE
  834. *
  835. ELSE IF( NOTRNA .AND. .NOT.NOTRNB ) THEN
  836. *
  837. * Solve A*X + ISGN*X*B**T = scale*C.
  838. *
  839. * The (K,L)th block of X is determined starting from
  840. * bottom-right corner column by column by
  841. *
  842. * A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L)**T = C(K,L) - R(K,L)
  843. *
  844. * Where
  845. * M N
  846. * R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)**T].
  847. * I=K+1 J=L+1
  848. *
  849. * Start column loop (index = L)
  850. * L1 (L2): column index of the first (last) row of X(K,L)
  851. *
  852. LNEXT = N
  853. DO 240 L = N, 1, -1
  854. IF( L.GT.LNEXT )
  855. $ GO TO 240
  856. IF( L.EQ.1 ) THEN
  857. L1 = L
  858. L2 = L
  859. ELSE
  860. IF( B( L, L-1 ).NE.ZERO ) THEN
  861. L1 = L - 1
  862. L2 = L
  863. LNEXT = L - 2
  864. ELSE
  865. L1 = L
  866. L2 = L
  867. LNEXT = L - 1
  868. END IF
  869. END IF
  870. *
  871. * Start row loop (index = K)
  872. * K1 (K2): row index of the first (last) row of X(K,L)
  873. *
  874. KNEXT = M
  875. DO 230 K = M, 1, -1
  876. IF( K.GT.KNEXT )
  877. $ GO TO 230
  878. IF( K.EQ.1 ) THEN
  879. K1 = K
  880. K2 = K
  881. ELSE
  882. IF( A( K, K-1 ).NE.ZERO ) THEN
  883. K1 = K - 1
  884. K2 = K
  885. KNEXT = K - 2
  886. ELSE
  887. K1 = K
  888. K2 = K
  889. KNEXT = K - 1
  890. END IF
  891. END IF
  892. *
  893. IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN
  894. SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
  895. $ C( MIN( K1+1, M ), L1 ), 1 )
  896. SUMR = DDOT( N-L1, C( K1, MIN( L1+1, N ) ), LDC,
  897. $ B( L1, MIN( L1+1, N ) ), LDB )
  898. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  899. SCALOC = ONE
  900. *
  901. A11 = A( K1, K1 ) + SGN*B( L1, L1 )
  902. DA11 = ABS( A11 )
  903. IF( DA11.LE.SMIN ) THEN
  904. A11 = SMIN
  905. DA11 = SMIN
  906. INFO = 1
  907. END IF
  908. DB = ABS( VEC( 1, 1 ) )
  909. IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN
  910. IF( DB.GT.BIGNUM*DA11 )
  911. $ SCALOC = ONE / DB
  912. END IF
  913. X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11
  914. *
  915. IF( SCALOC.NE.ONE ) THEN
  916. DO 190 J = 1, N
  917. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  918. 190 CONTINUE
  919. SCALE = SCALE*SCALOC
  920. END IF
  921. C( K1, L1 ) = X( 1, 1 )
  922. *
  923. ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN
  924. *
  925. SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
  926. $ C( MIN( K2+1, M ), L1 ), 1 )
  927. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  928. $ B( L1, MIN( L2+1, N ) ), LDB )
  929. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  930. *
  931. SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
  932. $ C( MIN( K2+1, M ), L1 ), 1 )
  933. SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
  934. $ B( L1, MIN( L2+1, N ) ), LDB )
  935. VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
  936. *
  937. CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, A( K1, K1 ),
  938. $ LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ),
  939. $ ZERO, X, 2, SCALOC, XNORM, IERR )
  940. IF( IERR.NE.0 )
  941. $ INFO = 1
  942. *
  943. IF( SCALOC.NE.ONE ) THEN
  944. DO 200 J = 1, N
  945. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  946. 200 CONTINUE
  947. SCALE = SCALE*SCALOC
  948. END IF
  949. C( K1, L1 ) = X( 1, 1 )
  950. C( K2, L1 ) = X( 2, 1 )
  951. *
  952. ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN
  953. *
  954. SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
  955. $ C( MIN( K1+1, M ), L1 ), 1 )
  956. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  957. $ B( L1, MIN( L2+1, N ) ), LDB )
  958. VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) )
  959. *
  960. SUML = DDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA,
  961. $ C( MIN( K1+1, M ), L2 ), 1 )
  962. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  963. $ B( L2, MIN( L2+1, N ) ), LDB )
  964. VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) )
  965. *
  966. CALL DLALN2( .FALSE., 2, 1, SMIN, ONE, B( L1, L1 ),
  967. $ LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ),
  968. $ ZERO, X, 2, SCALOC, XNORM, IERR )
  969. IF( IERR.NE.0 )
  970. $ INFO = 1
  971. *
  972. IF( SCALOC.NE.ONE ) THEN
  973. DO 210 J = 1, N
  974. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  975. 210 CONTINUE
  976. SCALE = SCALE*SCALOC
  977. END IF
  978. C( K1, L1 ) = X( 1, 1 )
  979. C( K1, L2 ) = X( 2, 1 )
  980. *
  981. ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN
  982. *
  983. SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
  984. $ C( MIN( K2+1, M ), L1 ), 1 )
  985. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  986. $ B( L1, MIN( L2+1, N ) ), LDB )
  987. VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR )
  988. *
  989. SUML = DDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA,
  990. $ C( MIN( K2+1, M ), L2 ), 1 )
  991. SUMR = DDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC,
  992. $ B( L2, MIN( L2+1, N ) ), LDB )
  993. VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR )
  994. *
  995. SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
  996. $ C( MIN( K2+1, M ), L1 ), 1 )
  997. SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
  998. $ B( L1, MIN( L2+1, N ) ), LDB )
  999. VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR )
  1000. *
  1001. SUML = DDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA,
  1002. $ C( MIN( K2+1, M ), L2 ), 1 )
  1003. SUMR = DDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC,
  1004. $ B( L2, MIN( L2+1, N ) ), LDB )
  1005. VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR )
  1006. *
  1007. CALL DLASY2( .FALSE., .TRUE., ISGN, 2, 2, A( K1, K1 ),
  1008. $ LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X,
  1009. $ 2, XNORM, IERR )
  1010. IF( IERR.NE.0 )
  1011. $ INFO = 1
  1012. *
  1013. IF( SCALOC.NE.ONE ) THEN
  1014. DO 220 J = 1, N
  1015. CALL DSCAL( M, SCALOC, C( 1, J ), 1 )
  1016. 220 CONTINUE
  1017. SCALE = SCALE*SCALOC
  1018. END IF
  1019. C( K1, L1 ) = X( 1, 1 )
  1020. C( K1, L2 ) = X( 1, 2 )
  1021. C( K2, L1 ) = X( 2, 1 )
  1022. C( K2, L2 ) = X( 2, 2 )
  1023. END IF
  1024. *
  1025. 230 CONTINUE
  1026. 240 CONTINUE
  1027. *
  1028. END IF
  1029. *
  1030. RETURN
  1031. *
  1032. * End of DTRSYL
  1033. *
  1034. END
  1035.  
  1036.  
  1037.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales