You can not select more than 25 topics Topics must start with a chinese character,a letter or number, can include dashes ('-') and can be up to 35 characters long.

clacn2.f 8.0 kB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298
  1. *> \brief \b CLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vector products.
  2. *
  3. * =========== DOCUMENTATION ===========
  4. *
  5. * Online html documentation available at
  6. * http://www.netlib.org/lapack/explore-html/
  7. *
  8. *> \htmlonly
  9. *> Download CLACN2 + dependencies
  10. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/clacn2.f">
  11. *> [TGZ]</a>
  12. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/clacn2.f">
  13. *> [ZIP]</a>
  14. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/clacn2.f">
  15. *> [TXT]</a>
  16. *> \endhtmlonly
  17. *
  18. * Definition:
  19. * ===========
  20. *
  21. * SUBROUTINE CLACN2( N, V, X, EST, KASE, ISAVE )
  22. *
  23. * .. Scalar Arguments ..
  24. * INTEGER KASE, N
  25. * REAL EST
  26. * ..
  27. * .. Array Arguments ..
  28. * INTEGER ISAVE( 3 )
  29. * COMPLEX V( * ), X( * )
  30. * ..
  31. *
  32. *
  33. *> \par Purpose:
  34. * =============
  35. *>
  36. *> \verbatim
  37. *>
  38. *> CLACN2 estimates the 1-norm of a square, complex matrix A.
  39. *> Reverse communication is used for evaluating matrix-vector products.
  40. *> \endverbatim
  41. *
  42. * Arguments:
  43. * ==========
  44. *
  45. *> \param[in] N
  46. *> \verbatim
  47. *> N is INTEGER
  48. *> The order of the matrix. N >= 1.
  49. *> \endverbatim
  50. *>
  51. *> \param[out] V
  52. *> \verbatim
  53. *> V is COMPLEX array, dimension (N)
  54. *> On the final return, V = A*W, where EST = norm(V)/norm(W)
  55. *> (W is not returned).
  56. *> \endverbatim
  57. *>
  58. *> \param[in,out] X
  59. *> \verbatim
  60. *> X is COMPLEX array, dimension (N)
  61. *> On an intermediate return, X should be overwritten by
  62. *> A * X, if KASE=1,
  63. *> A**H * X, if KASE=2,
  64. *> where A**H is the conjugate transpose of A, and CLACN2 must be
  65. *> re-called with all the other parameters unchanged.
  66. *> \endverbatim
  67. *>
  68. *> \param[in,out] EST
  69. *> \verbatim
  70. *> EST is REAL
  71. *> On entry with KASE = 1 or 2 and ISAVE(1) = 3, EST should be
  72. *> unchanged from the previous call to CLACN2.
  73. *> On exit, EST is an estimate (a lower bound) for norm(A).
  74. *> \endverbatim
  75. *>
  76. *> \param[in,out] KASE
  77. *> \verbatim
  78. *> KASE is INTEGER
  79. *> On the initial call to CLACN2, KASE should be 0.
  80. *> On an intermediate return, KASE will be 1 or 2, indicating
  81. *> whether X should be overwritten by A * X or A**H * X.
  82. *> On the final return from CLACN2, KASE will again be 0.
  83. *> \endverbatim
  84. *>
  85. *> \param[in,out] ISAVE
  86. *> \verbatim
  87. *> ISAVE is INTEGER array, dimension (3)
  88. *> ISAVE is used to save variables between calls to SLACN2
  89. *> \endverbatim
  90. *
  91. * Authors:
  92. * ========
  93. *
  94. *> \author Univ. of Tennessee
  95. *> \author Univ. of California Berkeley
  96. *> \author Univ. of Colorado Denver
  97. *> \author NAG Ltd.
  98. *
  99. *> \date December 2016
  100. *
  101. *> \ingroup complexOTHERauxiliary
  102. *
  103. *> \par Further Details:
  104. * =====================
  105. *>
  106. *> \verbatim
  107. *>
  108. *> Originally named CONEST, dated March 16, 1988.
  109. *>
  110. *> Last modified: April, 1999
  111. *>
  112. *> This is a thread safe version of CLACON, which uses the array ISAVE
  113. *> in place of a SAVE statement, as follows:
  114. *>
  115. *> CLACON CLACN2
  116. *> JUMP ISAVE(1)
  117. *> J ISAVE(2)
  118. *> ITER ISAVE(3)
  119. *> \endverbatim
  120. *
  121. *> \par Contributors:
  122. * ==================
  123. *>
  124. *> Nick Higham, University of Manchester
  125. *
  126. *> \par References:
  127. * ================
  128. *>
  129. *> N.J. Higham, "FORTRAN codes for estimating the one-norm of
  130. *> a real or complex matrix, with applications to condition estimation",
  131. *> ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988.
  132. *>
  133. * =====================================================================
  134. SUBROUTINE CLACN2( N, V, X, EST, KASE, ISAVE )
  135. *
  136. * -- LAPACK auxiliary routine (version 3.7.0) --
  137. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  138. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  139. * December 2016
  140. *
  141. * .. Scalar Arguments ..
  142. INTEGER KASE, N
  143. REAL EST
  144. * ..
  145. * .. Array Arguments ..
  146. INTEGER ISAVE( 3 )
  147. COMPLEX V( * ), X( * )
  148. * ..
  149. *
  150. * =====================================================================
  151. *
  152. * .. Parameters ..
  153. INTEGER ITMAX
  154. PARAMETER ( ITMAX = 5 )
  155. REAL ONE, TWO
  156. PARAMETER ( ONE = 1.0E0, TWO = 2.0E0 )
  157. COMPLEX CZERO, CONE
  158. PARAMETER ( CZERO = ( 0.0E0, 0.0E0 ),
  159. $ CONE = ( 1.0E0, 0.0E0 ) )
  160. * ..
  161. * .. Local Scalars ..
  162. INTEGER I, JLAST
  163. REAL ABSXI, ALTSGN, ESTOLD, SAFMIN, TEMP
  164. * ..
  165. * .. External Functions ..
  166. INTEGER ICMAX1
  167. REAL SCSUM1, SLAMCH
  168. EXTERNAL ICMAX1, SCSUM1, SLAMCH
  169. * ..
  170. * .. External Subroutines ..
  171. EXTERNAL CCOPY
  172. * ..
  173. * .. Intrinsic Functions ..
  174. INTRINSIC ABS, AIMAG, CMPLX, REAL
  175. * ..
  176. * .. Executable Statements ..
  177. *
  178. SAFMIN = SLAMCH( 'Safe minimum' )
  179. IF( KASE.EQ.0 ) THEN
  180. DO 10 I = 1, N
  181. X( I ) = CMPLX( ONE / REAL( N ) )
  182. 10 CONTINUE
  183. KASE = 1
  184. ISAVE( 1 ) = 1
  185. RETURN
  186. END IF
  187. *
  188. GO TO ( 20, 40, 70, 90, 120 )ISAVE( 1 )
  189. *
  190. * ................ ENTRY (ISAVE( 1 ) = 1)
  191. * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY A*X.
  192. *
  193. 20 CONTINUE
  194. IF( N.EQ.1 ) THEN
  195. V( 1 ) = X( 1 )
  196. EST = ABS( V( 1 ) )
  197. * ... QUIT
  198. GO TO 130
  199. END IF
  200. EST = SCSUM1( N, X, 1 )
  201. *
  202. DO 30 I = 1, N
  203. ABSXI = ABS( X( I ) )
  204. IF( ABSXI.GT.SAFMIN ) THEN
  205. X( I ) = CMPLX( REAL( X( I ) ) / ABSXI,
  206. $ AIMAG( X( I ) ) / ABSXI )
  207. ELSE
  208. X( I ) = CONE
  209. END IF
  210. 30 CONTINUE
  211. KASE = 2
  212. ISAVE( 1 ) = 2
  213. RETURN
  214. *
  215. * ................ ENTRY (ISAVE( 1 ) = 2)
  216. * FIRST ITERATION. X HAS BEEN OVERWRITTEN BY CTRANS(A)*X.
  217. *
  218. 40 CONTINUE
  219. ISAVE( 2 ) = ICMAX1( N, X, 1 )
  220. ISAVE( 3 ) = 2
  221. *
  222. * MAIN LOOP - ITERATIONS 2,3,...,ITMAX.
  223. *
  224. 50 CONTINUE
  225. DO 60 I = 1, N
  226. X( I ) = CZERO
  227. 60 CONTINUE
  228. X( ISAVE( 2 ) ) = CONE
  229. KASE = 1
  230. ISAVE( 1 ) = 3
  231. RETURN
  232. *
  233. * ................ ENTRY (ISAVE( 1 ) = 3)
  234. * X HAS BEEN OVERWRITTEN BY A*X.
  235. *
  236. 70 CONTINUE
  237. CALL CCOPY( N, X, 1, V, 1 )
  238. ESTOLD = EST
  239. EST = SCSUM1( N, V, 1 )
  240. *
  241. * TEST FOR CYCLING.
  242. IF( EST.LE.ESTOLD )
  243. $ GO TO 100
  244. *
  245. DO 80 I = 1, N
  246. ABSXI = ABS( X( I ) )
  247. IF( ABSXI.GT.SAFMIN ) THEN
  248. X( I ) = CMPLX( REAL( X( I ) ) / ABSXI,
  249. $ AIMAG( X( I ) ) / ABSXI )
  250. ELSE
  251. X( I ) = CONE
  252. END IF
  253. 80 CONTINUE
  254. KASE = 2
  255. ISAVE( 1 ) = 4
  256. RETURN
  257. *
  258. * ................ ENTRY (ISAVE( 1 ) = 4)
  259. * X HAS BEEN OVERWRITTEN BY CTRANS(A)*X.
  260. *
  261. 90 CONTINUE
  262. JLAST = ISAVE( 2 )
  263. ISAVE( 2 ) = ICMAX1( N, X, 1 )
  264. IF( ( ABS( X( JLAST ) ).NE.ABS( X( ISAVE( 2 ) ) ) ) .AND.
  265. $ ( ISAVE( 3 ).LT.ITMAX ) ) THEN
  266. ISAVE( 3 ) = ISAVE( 3 ) + 1
  267. GO TO 50
  268. END IF
  269. *
  270. * ITERATION COMPLETE. FINAL STAGE.
  271. *
  272. 100 CONTINUE
  273. ALTSGN = ONE
  274. DO 110 I = 1, N
  275. X( I ) = CMPLX( ALTSGN*( ONE + REAL( I-1 ) / REAL( N-1 ) ) )
  276. ALTSGN = -ALTSGN
  277. 110 CONTINUE
  278. KASE = 1
  279. ISAVE( 1 ) = 5
  280. RETURN
  281. *
  282. * ................ ENTRY (ISAVE( 1 ) = 5)
  283. * X HAS BEEN OVERWRITTEN BY A*X.
  284. *
  285. 120 CONTINUE
  286. TEMP = TWO*( SCSUM1( N, X, 1 ) / REAL( 3*N ) )
  287. IF( TEMP.GT.EST ) THEN
  288. CALL CCOPY( N, X, 1, V, 1 )
  289. EST = TEMP
  290. END IF
  291. *
  292. 130 CONTINUE
  293. KASE = 0
  294. RETURN
  295. *
  296. * End of CLACN2
  297. *
  298. END