1c----------------------------------------------------------------------- 2c\BeginDoc 3c 4c\Name: sngets 5c 6c\Description: 7c Given the eigenvalues of the upper Hessenberg matrix H, 8c computes the NP shifts AMU that are zeros of the polynomial of 9c degree NP which filters out components of the unwanted eigenvectors 10c corresponding to the AMU's based on some given criteria. 11c 12c NOTE: call this even in the case of user specified shifts in order 13c to sort the eigenvalues, and error bounds of H for later use. 14c 15c\Usage: 16c call sngets 17c ( ISHIFT, WHICH, KEV, NP, RITZR, RITZI, BOUNDS, SHIFTR, SHIFTI ) 18c 19c\Arguments 20c ISHIFT Integer. (INPUT) 21c Method for selecting the implicit shifts at each iteration. 22c ISHIFT = 0: user specified shifts 23c ISHIFT = 1: exact shift with respect to the matrix H. 24c 25c WHICH Character*2. (INPUT) 26c Shift selection criteria. 27c 'LM' -> want the KEV eigenvalues of largest magnitude. 28c 'SM' -> want the KEV eigenvalues of smallest magnitude. 29c 'LR' -> want the KEV eigenvalues of largest real part. 30c 'SR' -> want the KEV eigenvalues of smallest real part. 31c 'LI' -> want the KEV eigenvalues of largest imaginary part. 32c 'SI' -> want the KEV eigenvalues of smallest imaginary part. 33c 34c KEV Integer. (INPUT/OUTPUT) 35c INPUT: KEV+NP is the size of the matrix H. 36c OUTPUT: Possibly increases KEV by one to keep complex conjugate 37c pairs together. 38c 39c NP Integer. (INPUT/OUTPUT) 40c Number of implicit shifts to be computed. 41c OUTPUT: Possibly decreases NP by one to keep complex conjugate 42c pairs together. 43c 44c RITZR, Real array of length KEV+NP. (INPUT/OUTPUT) 45c RITZI On INPUT, RITZR and RITZI contain the real and imaginary 46c parts of the eigenvalues of H. 47c On OUTPUT, RITZR and RITZI are sorted so that the unwanted 48c eigenvalues are in the first NP locations and the wanted 49c portion is in the last KEV locations. When exact shifts are 50c selected, the unwanted part corresponds to the shifts to 51c be applied. Also, if ISHIFT .eq. 1, the unwanted eigenvalues 52c are further sorted so that the ones with largest Ritz values 53c are first. 54c 55c BOUNDS Real array of length KEV+NP. (INPUT/OUTPUT) 56c Error bounds corresponding to the ordering in RITZ. 57c 58c SHIFTR, SHIFTI *** USE deprecated as of version 2.1. *** 59c 60c 61c\EndDoc 62c 63c----------------------------------------------------------------------- 64c 65c\BeginLib 66c 67c\Local variables: 68c xxxxxx real 69c 70c\Routines called: 71c ssortc ARPACK sorting routine. 72c scopy Level 1 BLAS that copies one vector to another . 73c 74c\Author 75c Danny Sorensen Phuong Vu 76c Richard Lehoucq CRPC / Rice University 77c Dept. of Computational & Houston, Texas 78c Applied Mathematics 79c Rice University 80c Houston, Texas 81c 82c\Revision history: 83c xx/xx/92: Version ' 2.1' 84c 85c\SCCS Information: @(#) 86c FILE: ngets.F SID: 2.3 DATE OF SID: 4/20/96 RELEASE: 2 87c 88c\Remarks 89c 1. xxxx 90c 91c\EndLib 92c 93c----------------------------------------------------------------------- 94c 95 subroutine sngets ( ishift, which, kev, np, ritzr, ritzi, bounds, 96 & shiftr, shifti ) 97c 98c %----------------------------------------------------% 99c | Include files for debugging and timing information | 100c %----------------------------------------------------% 101c 102 include 'debug.h' 103 include 'stat.h' 104c 105c %------------------% 106c | Scalar Arguments | 107c %------------------% 108c 109 character*2 which 110 integer ishift, kev, np 111c 112c %-----------------% 113c | Array Arguments | 114c %-----------------% 115c 116 Real 117 & bounds(kev+np), ritzr(kev+np), ritzi(kev+np), 118 & shiftr(1), shifti(1) 119c 120c %------------% 121c | Parameters | 122c %------------% 123c 124 Real 125 & one, zero 126 parameter (one = 1.0, zero = 0.0) 127c 128c %---------------% 129c | Local Scalars | 130c %---------------% 131c 132 integer msglvl 133c 134c %----------------------% 135c | External Subroutines | 136c %----------------------% 137c 138 external scopy, ssortc, second 139c 140c %----------------------% 141c | Intrinsics Functions | 142c %----------------------% 143c 144 intrinsic abs 145c 146c %-----------------------% 147c | Executable Statements | 148c %-----------------------% 149c 150c %-------------------------------% 151c | Initialize timing statistics | 152c | & message level for debugging | 153c %-------------------------------% 154c 155 call second (t0) 156 msglvl = mngets 157c 158c %----------------------------------------------------% 159c | LM, SM, LR, SR, LI, SI case. | 160c | Sort the eigenvalues of H into the desired order | 161c | and apply the resulting order to BOUNDS. | 162c | The eigenvalues are sorted so that the wanted part | 163c | are always in the last KEV locations. | 164c | We first do a pre-processing sort in order to keep | 165c | complex conjugate pairs together | 166c %----------------------------------------------------% 167c 168 if (which .eq. 'LM') then 169 call ssortc ('LR', .true., kev+np, ritzr, ritzi, bounds) 170 else if (which .eq. 'SM') then 171 call ssortc ('SR', .true., kev+np, ritzr, ritzi, bounds) 172 else if (which .eq. 'LR') then 173 call ssortc ('LM', .true., kev+np, ritzr, ritzi, bounds) 174 else if (which .eq. 'SR') then 175 call ssortc ('SM', .true., kev+np, ritzr, ritzi, bounds) 176 else if (which .eq. 'LI') then 177 call ssortc ('LM', .true., kev+np, ritzr, ritzi, bounds) 178 else if (which .eq. 'SI') then 179 call ssortc ('SM', .true., kev+np, ritzr, ritzi, bounds) 180 end if 181c 182 call ssortc (which, .true., kev+np, ritzr, ritzi, bounds) 183c 184c %-------------------------------------------------------% 185c | Increase KEV by one if the ( ritzr(np),ritzi(np) ) | 186c | = ( ritzr(np+1),-ritzi(np+1) ) and ritz(np) .ne. zero | 187c | Accordingly decrease NP by one. In other words keep | 188c | complex conjugate pairs together. | 189c %-------------------------------------------------------% 190c 191 if ( ( ritzr(np+1) - ritzr(np) ) .eq. zero 192 & .and. ( ritzi(np+1) + ritzi(np) ) .eq. zero ) then 193 np = np - 1 194 kev = kev + 1 195 end if 196c 197 if ( ishift .eq. 1 ) then 198c 199c %-------------------------------------------------------% 200c | Sort the unwanted Ritz values used as shifts so that | 201c | the ones with largest Ritz estimates are first | 202c | This will tend to minimize the effects of the | 203c | forward instability of the iteration when they shifts | 204c | are applied in subroutine snapps. | 205c | Be careful and use 'SR' since we want to sort BOUNDS! | 206c %-------------------------------------------------------% 207c 208 call ssortc ( 'SR', .true., np, bounds, ritzr, ritzi ) 209 end if 210c 211 call second (t1) 212 tngets = tngets + (t1 - t0) 213c 214 if (msglvl .gt. 0) then 215 call ivout (logfil, 1, kev, ndigit, '_ngets: KEV is') 216 call ivout (logfil, 1, np, ndigit, '_ngets: NP is') 217 call svout (logfil, kev+np, ritzr, ndigit, 218 & '_ngets: Eigenvalues of current H matrix -- real part') 219 call svout (logfil, kev+np, ritzi, ndigit, 220 & '_ngets: Eigenvalues of current H matrix -- imag part') 221 call svout (logfil, kev+np, bounds, ndigit, 222 & '_ngets: Ritz estimates of the current KEV+NP Ritz values') 223 end if 224c 225 return 226c 227c %---------------% 228c | End of sngets | 229c %---------------% 230c 231 end 232