1 subroutine ortran(nm,n,low,igh,a,ort,z) 2c 3 integer i,j,n,kl,mm,mp,nm,igh,low,mp1 4 double precision a(nm,igh),ort(igh),z(nm,n) 5 double precision g 6c 7c this subroutine is a translation of the algol procedure ortrans, 8c num. math. 16, 181-204(1970) by peters and wilkinson. 9c handbook for auto. comp., vol.ii-linear algebra, 372-395(1971). 10c 11c this subroutine accumulates the orthogonal similarity 12c transformations used in the reduction of a real general 13c matrix to upper hessenberg form by orthes. 14c 15c on input 16c 17c nm must be set to the row dimension of two-dimensional 18c array parameters as declared in the calling program 19c dimension statement. 20c 21c n is the order of the matrix. 22c 23c low and igh are integers determined by the balancing 24c subroutine balanc. if balanc has not been used, 25c set low=1, igh=n. 26c 27c a contains information about the orthogonal trans- 28c formations used in the reduction by orthes 29c in its strict lower triangle. 30c 31c ort contains further information about the trans- 32c formations used in the reduction by orthes. 33c only elements low through igh are used. 34c 35c on output 36c 37c z contains the transformation matrix produced in the 38c reduction by orthes. 39c 40c ort has been altered. 41c 42c questions and comments should be directed to burton s. garbow, 43c mathematics and computer science div, argonne national laboratory 44c 45c this version dated august 1983. 46c 47c ------------------------------------------------------------------ 48c 49c .......... initialize z to identity matrix .......... 50 do 80 j = 1, n 51c 52 do 60 i = 1, n 53 60 z(i,j) = 0.0d0 54c 55 z(j,j) = 1.0d0 56 80 continue 57c 58 kl = igh - low - 1 59 if (kl .lt. 1) go to 200 60c .......... for mp=igh-1 step -1 until low+1 do -- .......... 61 do 140 mm = 1, kl 62 mp = igh - mm 63 if (a(mp,mp-1) .eq. 0.0d0) go to 140 64 mp1 = mp + 1 65c 66 do 100 i = mp1, igh 67 100 ort(i) = a(i,mp-1) 68c 69 do 130 j = mp, igh 70 g = 0.0d0 71c 72 do 110 i = mp, igh 73 110 g = g + ort(i) * z(i,j) 74c .......... divisor below is negative of h formed in orthes. 75c double division avoids possible underflow .......... 76 g = (g / ort(mp)) / a(mp,mp-1) 77c 78 do 120 i = mp, igh 79 120 z(i,j) = z(i,j) + g * ort(i) 80c 81 130 continue 82c 83 140 continue 84c 85 200 return 86 end 87