given an m by n matrix a, this subroutine computes aq where q is the product of 2(n - 1) transformations
gv(n-1)*...*gv(1)*gw(1)*...*gw(n-1)
and gv(i), gw(i) are givens rotations in the (i,n) plane which eliminate elements in the i-th and n-th planes, respectively. q itself is not given, rather the information to recover the gv, gw rotations is supplied.
the subroutine statement is
subroutine r1mpyq(m,n,a,lda,v,w)
where
m is a positive integer input variable set to the number of rows of a.
n is a positive integer input variable set to the number of columns of a.
a is an m by n array. on input a must contain the matrix to be postmultiplied by the orthogonal matrix q described above. on output a*q has replaced a.
lda is a positive integer input variable not less than m which specifies the leading dimension of the array a.
v is an input array of length n. v(i) must contain the information necessary to recover the givens rotation gv(i) described above.
w is an input array of length n. w(i) must contain the information necessary to recover the givens rotation gw(i) described above.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer | :: | m | ||||
integer | :: | n | ||||
real(kind=wp) | :: | a(lda,n) | ||||
integer | :: | lda | ||||
real(kind=wp) | :: | v(n) | ||||
real(kind=wp) | :: | w(n) |
subroutine r1mpyq(m,n,a,lda,v,w) implicit none integer m , n , lda real(wp) a(lda,n) , v(n) , w(n) integer i , j , nmj , nm1 real(wp) cos , sin , temp ! ! apply the first set of givens rotations to a. ! nm1 = n - 1 if ( nm1>=1 ) then do nmj = 1 , nm1 j = n - nmj if ( abs(v(j))>one ) cos = one/v(j) if ( abs(v(j))>one ) sin = sqrt(one-cos**2) if ( abs(v(j))<=one ) sin = v(j) if ( abs(v(j))<=one ) cos = sqrt(one-sin**2) do i = 1 , m temp = cos*a(i,j) - sin*a(i,n) a(i,n) = sin*a(i,j) + cos*a(i,n) a(i,j) = temp enddo enddo ! ! apply the second set of givens rotations to a. ! do j = 1 , nm1 if ( abs(w(j))>one ) cos = one/w(j) if ( abs(w(j))>one ) sin = sqrt(one-cos**2) if ( abs(w(j))<=one ) sin = w(j) if ( abs(w(j))<=one ) cos = sqrt(one-sin**2) do i = 1 , m temp = cos*a(i,j) + sin*a(i,n) a(i,n) = -sin*a(i,j) + cos*a(i,n) a(i,j) = temp enddo enddo endif end subroutine r1mpyq