r1mpyq Subroutine

public subroutine r1mpyq(m, n, a, lda, v, w)

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.

Arguments

Type IntentOptional Attributes Name
integer :: m
integer :: n
real(kind=wp) :: a(lda,n)
integer :: lda
real(kind=wp) :: v(n)
real(kind=wp) :: w(n)

Called by

proc~~r1mpyq~~CalledByGraph proc~r1mpyq minpack_module::r1mpyq proc~hybrd minpack_module::hybrd proc~hybrd->proc~r1mpyq proc~hybrj minpack_module::hybrj proc~hybrj->proc~r1mpyq proc~hybrd1 minpack_module::hybrd1 proc~hybrd1->proc~hybrd proc~hybrj1 minpack_module::hybrj1 proc~hybrj1->proc~hybrj proc~halo_to_rv_diffcorr halo_orbit_module::halo_to_rv_diffcorr proc~halo_to_rv_diffcorr->proc~hybrd1

Source Code

    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