srlcal.f 4.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115
  1. *DECK SRLCAL
  2. SUBROUTINE SRLCAL (N, KMP, LL, MAXL, V, Q, RL, SNORMW, PROD,
  3. + R0NRM)
  4. C***BEGIN PROLOGUE SRLCAL
  5. C***SUBSIDIARY
  6. C***PURPOSE Internal routine for SGMRES.
  7. C***LIBRARY SLATEC (SLAP)
  8. C***CATEGORY D2A4, D2B4
  9. C***TYPE SINGLE PRECISION (SRLCAL-S, DRLCAL-D)
  10. C***KEYWORDS GENERALIZED MINIMUM RESIDUAL, ITERATIVE PRECONDITION,
  11. C NON-SYMMETRIC LINEAR SYSTEM, SLAP, SPARSE
  12. C***AUTHOR Brown, Peter, (LLNL), pnbrown@llnl.gov
  13. C Hindmarsh, Alan, (LLNL), alanh@llnl.gov
  14. C Seager, Mark K., (LLNL), seager@llnl.gov
  15. C Lawrence Livermore National Laboratory
  16. C PO Box 808, L-60
  17. C Livermore, CA 94550 (510) 423-3141
  18. C***DESCRIPTION
  19. C This routine calculates the scaled residual RL from the
  20. C V(I)'s.
  21. C *Usage:
  22. C INTEGER N, KMP, LL, MAXL
  23. C REAL V(N,LL), Q(2*MAXL), RL(N), SNORMW, PROD, R0NORM
  24. C
  25. C CALL SRLCAL(N, KMP, LL, MAXL, V, Q, RL, SNORMW, PROD, R0NRM)
  26. C
  27. C *Arguments:
  28. C N :IN Integer
  29. C The order of the matrix A, and the lengths
  30. C of the vectors SR, SZ, R0 and Z.
  31. C KMP :IN Integer
  32. C The number of previous V vectors the new vector VNEW
  33. C must be made orthogonal to. (KMP .le. MAXL)
  34. C LL :IN Integer
  35. C The current dimension of the Krylov subspace.
  36. C MAXL :IN Integer
  37. C The maximum dimension of the Krylov subspace.
  38. C V :IN Real V(N,LL)
  39. C The N x LL array containing the orthogonal vectors
  40. C V(*,1) to V(*,LL).
  41. C Q :IN Real Q(2*MAXL)
  42. C A real array of length 2*MAXL containing the components
  43. C of the Givens rotations used in the QR decomposition
  44. C of HES. It is loaded in SHEQR and used in SHELS.
  45. C RL :OUT Real RL(N)
  46. C The residual vector RL. This is either SB*(B-A*XL) if
  47. C not preconditioning or preconditioning on the right,
  48. C or SB*(M-inverse)*(B-A*XL) if preconditioning on the
  49. C left.
  50. C SNORMW :IN Real
  51. C Scale factor.
  52. C PROD :IN Real
  53. C The product s1*s2*...*sl = the product of the sines of the
  54. C Givens rotations used in the QR factorization of
  55. C the Hessenberg matrix HES.
  56. C R0NRM :IN Real
  57. C The scaled norm of initial residual R0.
  58. C
  59. C***SEE ALSO SGMRES
  60. C***ROUTINES CALLED SCOPY, SSCAL
  61. C***REVISION HISTORY (YYMMDD)
  62. C 871001 DATE WRITTEN
  63. C 881213 Previous REVISION DATE
  64. C 890915 Made changes requested at July 1989 CML Meeting. (MKS)
  65. C 890922 Numerous changes to prologue to make closer to SLATEC
  66. C standard. (FNF)
  67. C 890929 Numerous changes to reduce SP/DP differences. (FNF)
  68. C 910411 Prologue converted to Version 4.0 format. (BAB)
  69. C 910506 Made subsidiary to SGMRES. (FNF)
  70. C 920511 Added complete declaration section. (WRB)
  71. C***END PROLOGUE SRLCAL
  72. C The following is for optimized compilation on LLNL/LTSS Crays.
  73. CLLL. OPTIMIZE
  74. C .. Scalar Arguments ..
  75. REAL PROD, R0NRM, SNORMW
  76. INTEGER KMP, LL, MAXL, N
  77. C .. Array Arguments ..
  78. REAL Q(*), RL(N), V(N,*)
  79. C .. Local Scalars ..
  80. REAL C, S, TEM
  81. INTEGER I, I2, IP1, K, LLM1, LLP1
  82. C .. External Subroutines ..
  83. EXTERNAL SCOPY, SSCAL
  84. C***FIRST EXECUTABLE STATEMENT SRLCAL
  85. IF (KMP .EQ. MAXL) THEN
  86. C
  87. C calculate RL. Start by copying V(*,1) into RL.
  88. C
  89. CALL SCOPY(N, V(1,1), 1, RL, 1)
  90. LLM1 = LL - 1
  91. DO 20 I = 1,LLM1
  92. IP1 = I + 1
  93. I2 = I*2
  94. S = Q(I2)
  95. C = Q(I2-1)
  96. DO 10 K = 1,N
  97. RL(K) = S*RL(K) + C*V(K,IP1)
  98. 10 CONTINUE
  99. 20 CONTINUE
  100. S = Q(2*LL)
  101. C = Q(2*LL-1)/SNORMW
  102. LLP1 = LL + 1
  103. DO 30 K = 1,N
  104. RL(K) = S*RL(K) + C*V(K,LLP1)
  105. 30 CONTINUE
  106. ENDIF
  107. C
  108. C When KMP < MAXL, RL vector already partially calculated.
  109. C Scale RL by R0NRM*PROD to obtain the residual RL.
  110. C
  111. TEM = R0NRM*PROD
  112. CALL SSCAL(N, TEM, RL, 1)
  113. RETURN
  114. C------------- LAST LINE OF SRLCAL FOLLOWS ----------------------------
  115. END