formt Subroutine

private subroutine formt(m, Wt, Sy, Ss, Col, Theta, Info)

This subroutine forms the upper half of the pos. def. and symm. T = thetaSS + LD^(-1)L', stores T in the upper triangle of the array wt, and performs the Cholesky factorization of T to produce JJ', with J' stored in the upper triangle of wt.

Credits

  • NEOS, November 1994. (Latest revision June 1996.) Optimization Technology Center. Argonne National Laboratory and Northwestern University. Written by Ciyou Zhu in collaboration with R.H. Byrd, P. Lu-Chen and J. Nocedal.

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: m
real(kind=wp) :: Wt(m,m)
real(kind=wp) :: Sy(m,m)
real(kind=wp) :: Ss(m,m)
integer :: Col
real(kind=wp), intent(in) :: Theta
integer :: Info

Calls

proc~~formt~~CallsGraph proc~formt lbfgsb_module::formt proc~dpofa lbfgsb_linpack_module::dpofa proc~formt->proc~dpofa proc~ddot lbfgsb_blas_module::ddot proc~dpofa->proc~ddot

Called by

proc~~formt~~CalledByGraph proc~formt lbfgsb_module::formt proc~mainlb lbfgsb_module::mainlb proc~mainlb->proc~formt proc~setulb lbfgsb_module::setulb proc~setulb->proc~mainlb

Source Code

      subroutine formt(m,Wt,Sy,Ss,Col,Theta,Info)
      implicit none

      integer,intent(in) :: m
      integer :: Col
      integer :: Info
      real(wp),intent(in) :: Theta
      real(wp) :: Wt(m,m)
      real(wp) :: Sy(m,m)
      real(wp) :: Ss(m,m)

      integer :: i , j , k , k1
      real(wp) :: ddum

      ! Form the upper half of  T = theta*SS + L*D^(-1)*L',
      ! store T in the upper triangle of the array wt.

      do j = 1 , Col
         Wt(1,j) = Theta*Ss(1,j)
      enddo
      do i = 2 , Col
         do j = i , Col
            k1 = min(i,j) - 1
            ddum = zero
            do k = 1 , k1
               ddum = ddum + Sy(i,k)*Sy(j,k)/Sy(k,k)
            enddo
            Wt(i,j) = ddum + Theta*Ss(i,j)
         enddo
      enddo

      ! Cholesky factorize T to J*J' with
      ! J' stored in the upper triangle of wt.

      call dpofa(Wt,m,Col,Info)
      if ( Info/=0 ) Info = -3

      end subroutine formt